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Reservoir operation plays an important role in managing a regional water 
distribution system. During normal periods of operation, when inflows are plentiful, it is 
not difficult to distribute the available stored water from different reservoirs to meet the 
planned demands imposed by competing users. However, during periods of drought, or 
when anticipating a drought, the planned demands cannot be met completely, and a water 
shortage occurs. Therefore, a standard policy is needed to provide the rules for operating 
the reservoir to efficiently utilize the water. As a result, the rule curves and hedging rules 
are introduced as guidelines for the reservoir management and operation. 
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This study develops a mixed-integer linear programming (MILP) model that 
simultaneously considers both the rule curves and hedging rules to optimize the 
management and operation for a multi-purpose, multi-reservoir system. To minimize the 
impact of drought, the hedging rules effectively reduce the ongoing water supply to 
balance with the target storage requirement. 

A mixed-integer nonlinear programming (MINLP) model is then developed based 
on the proposed MILP model to evaluate and modify the current rule curves and hedging 
rules. Due to the mathematical complexity of this MINLP model, a number of 
methodologies are proposed to solve this model. These methodologies are: Generalized 
Benders Decomposition, Penalty Coefficient Methods, Pseudo-integer Method, 
Simulated Annealing, and one mathematical transformation technique. The results of 
these methodologies are explored and investigated. After examining the performance of 
these methodologies, the transformation technique seems best suited to solving this 
proposed MINLP model. 

A multi-reservoir system in the southern region of Taiwan is used as an actual 
case study to demonstrate the capability and utility of the proposed models. The results 
obtained indicate that the MILP model can efficiently utilize the water resources when 
incorporating the current rule curves and hedging rules. In addition, the new set of rule 
curves and hedging rules obtained by the MINLP model can successfully improve the 
performance of the current system operation. 
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CHAPTER 1 INTRODUCTION 



In water resources planning and management, the general purpose of system 
operation for a regional water distribution system is to utilize water appropriately and 
efficiently. A regional water distribution system includes: rivers, streams, reservoirs, 
power plants, diversion structures, pumping stations, water users, channels and pipes. 
System operation must consider all the above components in order for the system to 
operate as efficiently as possible. Consequently, system managers need to take into 
account the system complexities involved to maximize efficiency. In general, the primary 
concern in system operation is water supply. Among system components, reservoir 
operation plays a significant role in fulfilling water supply. 

Reservoirs are usually constructed for two purposes: water storage and water 
release. Unfortunately, these two purposes conflict with each other. The decision between 
these two competing purposes of either water storage or release can severely affect water 
usage for the entire system. If the current amount of water release is abundant, the water 
for future use may not be plentiful and have serious impact during periods of drought. If 
the current stored water is plentiful, the reservoir may not have enough space to store 
potential flood water and cause severe downstream damage. Therefore, system managers 
are responsible for determining a standard operating guideline to appropriately operate 
reservoirs. 
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A large-scale regional water supply system normally includes several single- or 
multi-purpose reservoirs. A single-purpose reservoir is designed for a single purpose such 
as municipal water supply. In contrast, a multi-purpose reservoir has multiple functions 
including: hydropower generation, flood control, navigation, fish and wildlife 
enhancement, and recreation. In addition, a multi-purpose reservoir functions to control 
agricultural, municipal and industrial water supply. The complexities of a multi-purpose, 
multi-reservoir system generally require release decisions made through the use of an 
assisting tool such as a simulation or optimization model. 

In general, simulation models are able to capture the actual physical 
characteristics of a system. However, the process of obtaining the most beneficial 
operating policies involves the repeated implementation of simulation models through a 
trial and error process. In contrast, while optimization models may not consider all the 
complex details of a system (i.e. assumptions may be set beforehand), they are more 
efficient in providing optimal policies for decision making. 

In recent years, optimization models have been used successfully to manage and 
operate complex reservoir systems. Yeh (1985), Simonovic (1992), Wurbs (1993), Wurbs 
(1996), ReVelle (1997), and Momoh et al. (1999a,b) have reviewed optimization model 
literature extensively and provide evaluations on various models. The choice of an 
appropriate optimization model depends on the following considerations: system 
characteristics, data availability, and the specified objectives and constraints. 
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1.1 Scope of Study 

During normal periods of operation when inflow is plentiful, meeting reservoir 
operating goals and target storage does not pose a problem. However, during periods of 
drought or when anticipating a drought, it may not be possible to reach the storage target 
because of an inflow shortage. To minimize the impact of drought and consequent 
shortages in current and future water supply, reservoir rule curves are coupled with 
hedging rules to balance the shortage in water supply with the storage target. 

This study is aimed at developing a multi-purpose and multi-resrvoir model to 
optimize the joint operation for a regional water supply system with regard to reservoir 
rule curves and hedging rules. This proposed optimization model can help plan for the 
operation of a multi-reservoir system during drought periods. Furthermore, this 
optimization model will also be expanded to evaluate and modify current rule curves and 
hedging rules for superior operation performance. 

1.2 Objectives of Study 

The objectives of this study are: 
1 . to develop a mixed-integer linear programming (MILP) model to optimize the joint 
operation of a multi-purpose and multi-reservoir system, incorporating rule curves 
and hedging rules; 
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2. to expand the above-mentioned MILP model to a mixed-integer nonlinear 
programming (MINLP) model to evaluate and modify the current rule curves and 
hedging rules; and 

3. to apply the developed optimization models to a regional water distribution system in 
Southern Taiwan. 

13 Dissertation Overview 

The remaining content of this dissertation is outlined as follows. Chapter 2 
includes a literature review concerning the optimization of regional water distribution 
systems. In addition, this chapter also discusses a number of past studies concerning 
reservoir operation in normal and dry situations. Chapter 3 describes the basic concept of 
developing an optimization model to deal with the operation of regional water 
distribution systems. The introduction of reservoir rule curves and hedging rules is then 
presented to provide the basic idea for operating reservoirs. 

Chapter 4 presents a mixed-integer linear programming (MILP) model that is 
capable of optimizing the reservoir operation incorporating rule curves and hedging rules. 
The capability of this model is tested by a simplified system application. Chapter 5 
proposes a mixed-integer nonlinear programming (MINLP) model, which is used to 
obtain a new set of rule curves and hedging rules to improve current system operation. 
Several methodologies are presented to solve this model and the results are compared and 
discussed. The model is then applied to a simplified system to test its capability. 
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Chapter 6 describes the application of the proposed MILP and MINLP models on 
a regional water distribution system - the southern regional water supply system in 
Taiwan. The model results are then examined and discussed. Finally, the conclusions and 
future work of this study are presented in Chapter 7. At the end of each chapter, 
corresponding tables and figures are included for the convenience of the reader. 
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CHAPTER 2 LITERATURE REVIEW 



There are a lot of studies concerning the optimization of regional water 
distribution systems in literature. The following section will review a number of these 
studies. In addition, the subsequent section will also summarize several studies involving 
reservoir operation during normal and dry periods. 

2.1 Optimization of Regional Water Distribution Systems 

Optimization techniques have already proven their applicability and robustness to 
many subjects concerning regional water distribution systems. The application of 
optimization techniques have been demonstrated in the following areas: (1) optimal 
network design (Ostfeld and Shamir, 1996; Taher and Labadie, 1996; Varma et al., 1997; 
Xu and Goulter, 1999); (2) reliability analysis (Su et al., 1987; Guercio and Zuxin, 1997); 
(3) water network leakage analysis (Vairavamoorthy and Lumbers, 1998; Alonso et al., 
2000; Bach et al., 2000). 

The following studies focused on water supply in a regional water distribution 
system. Kuczera (1989) presented a network linear program formulation for determining 
water assignments in a multireservoir system over some time horizon. The computing 
efficiency of this formulation was demonstrated as superior to the conventional linear 
programming formulation. Hsu and Yeh (1992) developed a linear programming (LP) 
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procedure utilizing the Dantzig- Wolfe decomposition method to solve a facility planning 
problem for a water supply system. Diba et al. (1995) presented a planning model 
incorporating a directed graph algorithm and an LP procedure to assist the system 
operation for a large-scale water distribution system. This model demonstrated its 
applicability and versatility via a real water distribution system. Sun et al. (1995) reported 
an embedded generalized network flow model representing a water distribution problem 
and discovered a solution via a fast network algorithm. Yang et al. (1996a) proposed a 
mechanical reliability methodology using deterministic analysis to describe the 
probability of source-demand connectivity. Yang et al. (1996b) presented a performance 
reliability methodology using stochastic analysis to evaluate the supply level that a 
system can provide under random component failures. 

Occasionally, water quality is considered simultaneously with water supply in 
order to optimize the operation of a regional water distribution system. Ostfeld and 
Shamir (1993 a) reported a steady state, nonlinear operation model for a multiquality 
network. As opposed to other studies on directional graph, the undirectional pipes were 
incorporated in their methodology. The steady state optimization model was further 
extended by Ostfeld and Shamir (1993 b) to incorporate the unsteady state solute transport 
in the water distribution system. Yang et al. (2000) presented preliminary research on the 
optimization of a regional water distribution system with blending requirements. In their 
study, a nonlinear multicommodity flow model was proposed to accommodate the 
blending requirements and perfect mixing conditions. 
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2.2 Reservoir Operation 

In a regional water distribution system, the reservoir operation is an essential key 
in optimizing system practice. If inflow is plentiful, the reservoir operation is not a 
serious issue. However, this ideal situation may not happen in every water resource 
system. If the water is not plentiful, system planners and managers need to decide how to 
operate the reservoirs to reduce drought damage and utilize water more efficiently. 
Therefore, the development of a useful reservoir operation guideline is critical. 

In many practical situations, operating rules (also referred to as operating policies) 
are established at the planning stage of the proposed reservoir. These rules then provide 
guidelines for reservoir releases to meet planned demands. The following literature 
review will report past studies on reservoir operation in normal and drought situations. 

2.2.1 Normal Situation 

Bhaskar and Whitlatch (1980) developed a backward dynamic program algorithm 
to obtain the optimal monthly releases for a single multi-purpose reservoir. Karamouz 
and Houck (1982) proposed an algorithm that consisted of a deterministic dynamic 
program, a regression analysis, and a simulation model to generate the annual and 
monthly operating rules for a single reservoir. Sand (1984) used stochastic dynamic 
programming to obtain optimal operating policies for the water supply of reservoirs in 
parallel. Using a genetic algorithm, Oliveira and Loucks (1997) derived optimal 
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operating policies for a multi-reservoir system. The above-mentioned studies have 
focused on reaching the reservoir release target as closely as possible. 

2.2.2 Drought Situation 

A small number of studies on hedging rules have been reported in literature. Shih 
and ReVelle (1994) developed a continuous linear hedging rule for a single water supply 
reservoir. They formulated a nonlinear mixed integer programming model that minimized 
the maximum deficit while considering a constant demand. After obtaining the optimal 
rule, they converted the continuous hedging rule into multiple discrete hedging rules 
which are more appropriate for practical applications. Shih and ReVelle (1995) 
developed a multiple hedging model using mixed integer programming for a single water 
supply reservoir. The primary purpose of the model was to maximize the number of 
months in which no water rationing is required. Neelakantan and Pundarikanthan (2000) 
presented a simulation-optimization methodology using neural network and multiple 
hedging rules to improve the reservoir operation performance for a drinking water 
reservoir system. Although the case study of their research comprised of several 
reservoirs, the storages of these reservoirs were lumped together into an equivalent 
reservoir in the mathematical model. All of the previous published papers have focused 
on operation of a single reservoir system, the operation for a multi-reservoir system was 
not considered. 
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CHAPTER 3 REGIONAL WATER SUPPLY 
OPTIMIZATION PROBLEM 



Developing an optimization model to resolve a water distribution system problem 
requires transforming the physical system into a mathematical model. A network 
representation is used as the transformation approach and helps develop an optimization 
model. Since the guideline for operating reservoirs is crucial in operating a regional water 
supply system, it is very important to consider the guideline in a system optimization 
problem. The preliminary introduction of the guideline - rule curves and hedging rules - 
is stated in the subsequent section. 

3.1 Network Representation 

To facilitate analysis of the interrelationship between various system components, 
a large-scale water supply distribution system can be represented as a network structure 
which is a node-link combination. For instance: inflows, reservoirs, diversion structures, 
and demands are represented as nodes; while rivers, channels and pipes are represented as 
arcs (Figure 3-1). A network consists of a set of nodes linked by arcs. The benefit of 
using a network to reconfigure the original water supply system is the ease with which 
the system can be clarified and formulated by mathematical models. 
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In general, the optimization of a regional water distribution system can be 
formulated as the following minimum-cost problem: 

r 

Min - Z Z c o.^ • x (..^ (3_1) 

'•'■ l x i,.u~ I *(„„=*„ • V/' eN (3-2) 

0,/)eA (/.,)eA 

x~ M Si w ^x- , V(/,y)eA (3-3) 

where 

/ = time index; 

T = total number of time periods; 

i,j = node index; 

(ij) = arc that emanates from node / and terminates at node/; 

X(ij).t - nonnegative flow in arc (ij); 

X o"jU = l0Wer b0Und 0f Xfij).h 

x uju = u PP er bound of x rw .,; 

c M./ = weighting factor (unit cost) for arc (7,/); 

b hl = source/sink term at node/; 

N = node set of the network; 

A = arc set of the network. 

In this formulation, the original system elements are characterized by nodes and 
arcs. Each node (or arc) can represent one specific component with a corresponding 
physical characteristic in the system. The decision variables in the optimization model are 
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the directional flows in each of the arcs in the configuration. In the composite objective 
function, Eq. (3-1), the weighting factor for each x^,, reflects the priority of each 
objective. In general, the determination of the weighting factor of each objective requires 
a systematic analysis of system characteristics to appropriately reflect the preference of 
the operator. Eq. (3-2) is the general form of a set of continuity equations which can be 
applied to inflow, diversion, junction, demand or reservoir nodes. Eq. (3-3) specifies the 
lower and upper bounds of each nonnegative flow variable. 

Figure 3-1 shows a network configuration at one time period. For planning 
purposes, the multi-time-period analysis needs to be taken into account. In the network 
representation, the storage of a reservoir is represented as an arc, termed the "reservoir 
carry-over arc", and this arc connects two adjacent time periods in a planning horizon 
(Figure 3-2) to construct a large-scale network structure for optimization. Continuity at 
the reservoir node is maintained by the continuity equation at the node. Any nodes with a 
carry-over storage function other than that of reservoirs can be represented similarly. In a 
large-scale network structure, each node or arc at each time period is considered 
simultaneously with corresponding data to implement system optimization. 

3.2 Reservoir Rule Curves and Hedging Rules 

Figure 3-3 shows the rule curves associated with a typical multi-purpose 
reservoir. In this example, the three curves are the flood control curve (the upper curve), 
the target storage curve (the middle curve), and the firm storage curve (the lower curve). 
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These curves are determined at the planning stage in order to provide guidelines for 
operating the reservoir. Although rule curves usually remain unchanged from year to 
year, they can, if necessary, be updated. In general, the flood control reservation is 
strictly observed. During normal periods of operation, all planned demands are met at 
100%, and the reservoir storage is kept at or above the established target storage level. 
However, during periods of drought when the reservoir storage falls below the target 
storage curve, the reservoir release is reduced to ensure that a sufficient amount of water 
remains in storage for future water supply. This situation creates a trade-off between 
meeting ongoing demands and maintaining adequate reservoir storage when inflow is 
insufficient. Using hedging rules along with the rule curves, a balance is achieved 
between water supply shortage and the storage objectives. 

A discussion on the use of rule curves and hedging rules is provided herein. 
Considering that during a specific time period, if the beginning reservoir storage is 
greater than that of the target storage curve (in zone 3), all planned demands are met at 
100%. If the beginning reservoir storage is greater than that of the firm storage curve and 
less than that of the target storage curve (in zone 2), the reservoir release for meeting the 
planned demand is reduced, for instance, by 25%; that is, the original demand is met at its 
75% level. If the beginning reservoir storage is less than that of the firm storage curve (in 
zone 1), then the reservoir release for meeting the planned demand is reduced even 
further, for example, by 50%; that is, the original demand is met at its 50% level. In other 
words, the lower the reservoir storage zone from which water is released to satisfy the 
planned demand, the higher the reduction of the demand. 
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Figure 3-1 An Example of Water Resources Network Configuration 
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Figure 3-3 Rule Curves for a Multi-purpose Reservoir 
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CHAPTER 4 MIXED-INTEGER LINEAR 
PROGRAMMING MODEL 



Sun et al. (1995) formulated a large-scale regional water supply problem in terms 
of an embedded generalized network flow model and discovered a solution through the 
use of a fast network algorithm. The model was for planning purposes, and hedging was 
not considered. To consider hedging, Yeh and Yang (1997) developed a network flow 
model in which the reservoir storage is divided into several volumes; these volumes are 
incorporated in the objective function as penalty terms and the associated cost 
coefficients are systematically determined by a sensitivity analysis. Hsu and Cheng 
(2002) used a similar approach and developed a network flow model for the Tanshui river 
basin in Taiwan. However, this approach is basically applicable to a single-time-period 
optimization model, and does not assure that the optimal solution of the multi-time- 
period model matches the rule curves and the hedging rules. 

In this study, a mixed-integer linear programming (MILP) formulation is 
developed to incorporate the rule curves and hedging rules in a multi-purpose and multi- 
time-period optimization model for a regional water supply system. 
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4.1 Single Reservoir System 

A single reservoir system is first explored with specified rule curves and hedging 
rules. The continuity equation for the reservoir can be written as the following: 

S,+I,-0, -E, =S, +/ (4-1) 

S mmJ < S, < S maxJ (4-2) 

where 

S, = beginning storage; 
S,+i = ending storage; 
/, = inflow during time period t; 
O, = release during time period t; 
E, = evaporation loss during time period t; 
Smmj = minimum storage; 
Smax.t = maximum storage. 

During each time period /, the relationship between the rule curves and the 
hedging rules is shown as a step function (Figure 4-1), which determines the preferred 
amount of water supply in the system. This definition is revised from Tu et al. (2002) in 
which the upper bound of a supply link is determined by hedging. The following 
equations represent the relationship: 
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If S mmJ <S,< S firmJ , then P,=a r D, (4-3) 

If S fimJ <S,< S^, , then P, = a 2 ■ D, (4-4) 

If S large ,, ^S,< S^, , then P, = D, (4-5) 

where 

Sfirm. t = firm storage; 
Stargeu= target storage; 
D, = planned demand; 
P, = water supply; 
a t and a.? = rationing factors, and < a/ < 02 < 1 . 

The rationing factors, aj and a?, are determined at the planning stage, but can be 
modified to meet operational purposes. 

In the network configuration, every demand node has a corresponding link that 
terminates at the demand node and the flow in the corresponding link is the water supply 
delivery to the demand node. Hence, the water supply variable P t , a decision variable in 
the optimization model, represents the amount of water that flows into a demand node. 

In the above formulation, only one of the constraints in Eqs. (4-3), (4-4), and (4-5) 
can be active, and the active constraint represents the actual reservoir storage and its 
corresponding hedging. In short, the correct solution should observe both the rule curves 
and the hedging rules. It is important to note that the water supply is assumed as equal to 
the planned demand or reduced planned demand, depending on which hedging is active. 
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Incorporating the reservoir rule curves and hedging rules does not cause difficulty 
in a simulation model (Jain et al. 1998). However, it is necessary to introduce new 
variables and constraints to extend Eqs. (3-1) ~(3-3) to consider hedging in an 
optimization model. To deal with multiple hedging rules in the optimization model, the 
following constraints are introduced: 

S, =S U +S 2J +5 3i , (4-6) 

^■S mUtJ <S u <X u -S firmJ (4-7) 

tTj-SjHuZSvZAu-Sw (4-8) 

K-s^u^Sx^K-s^., (4-9) 

/> l =(Va l +^-a 2 +'0-D, (4-10) 

K+K+K= l (4-1 la) 

\,.<*2,A,e{0,7} (4-1 lb) 

where 

S, = reservoir storage, which is divided into 5/,,, £?,,, and S 3 ,, (nonnegative variables); 
J-i.i, An, and ky t - binary (0-1) variables. 

Eq. (4-11) ensures that only one of the As can be active so that the reservoir 
storage matches the hedging rules. For example, if A/., is active (i.e., A/,, = 1 , A?, = 0, X 3j = 
0), then only Eq. (4-7) is active, and Eqs. (4-8) and (4-9) are inactive. Therefore, the 
optimal solution guarantees that, if the beginning storage in time period t is in zone 1, the 
amount of the water supply can be obtained from Eq. (4-10) as a t fraction of A- 
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4.2 Multi-reservoir System 

The formulation for the rule curves and hedging rules of a multiple reservoir 
system can be extended from the concept described above. Figure 4-2 shows a two- 
reservoir system (R 1 and R 2 with storages S 1 and S 2 ) that concurrently supplies water to 
one demand node (D). These two reservoirs have their own rule curves and hedging rules. 
The total amount of water that supplies demand node D is the summation of the waters 
from R 1 and R 2 . Assume that the percentage from each reservoir, PE 1 and PE 2 , is 
determined beforehand. For example, PE 1 = 0.3 and PE 2 = 0.7. In the network 
configuration, D can be separated into two nodes, D 1 and D 2 , with their corresponding 
arcs, A 1 and A 2 . The capacity of A 1 is equal to PE 1 fraction of A and the capacity of A 2 is 
equal to PE fraction of A. Accordingly, the original configuration, Figure 4-2(a), can be 
transformed into a new configuration, Figure 4-2(b). In general, the rule curves and 
hedging rules for inter-relations between reservoirs can be represented by the following 
set of constraints: 

S* + /* - O* - £* = 5* , , V* e K (4-12) 

S,* = S k u + Si + Si , V* € K (4-13) 

NK 

A., ■ SL, * St, < A, ■ S k firm . t , V* e K (4-15) 

4, -5* m , < 5*, < 4 -S^ u , Vk e K (4-16) 
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/>/ =(4 -a* +4 a 2 * +4)-P£* -A , V*eK (4-18) 

tia + Aj + A,, = 1 ,VieK (4-19a) 

^,^2,'^e{0,7} , VAreK (4- 19b) 

where 

£ = reservoir index; 

K = subset of reservoirs in which the rule curves and hedging rules of each reservoir are 

concurrently considered; 
NK = number of reservoirs in K; 

. NK 

PEr = percentage of water that is supplied by reservoir k, ^ PE k = 1 . 

*=i 

In the case of a two reservoir system shown in Figure 4-2, K = {R 1 , R 2 }, and NK 
= 2, and k = 1,2. The sum of the waters that supply D 1 and D 2 is the total amount of water 
that supplies D. 

43 The Overall MILP Model 

The summary of the steady-state, MILP model for a multi-reservoir system is 
described as follows: 
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«6yec/ to £ x ( ,,,, - J x (ljU = Z>„ , V; € N (4-21) 

(jjHA (,,j)eA 

x t..jL = x l.ju + x l.».< + x l,u ' v ft » e A « ' V -/' 6 N w (4-22) 

.va: 
*(..„, = £<,), .VyeN D/ (4-23) 

4, • -C * *,',„, * 4, • 3L, • V (''- J) e A « , V* e K„ , V/, e N DI , Vy e N„ (4-24) 

4 • Sj-u < xj„, < 4, • S,^,, , V(i, y) e A fl , V* e K„ , ty e N DI , Vy e N„ (4-25) 

4, • C«, * <„, * 4, • 5L, , V(i, y) eA B , V* e K„ , V|i e N 0/ , Vy e N w (4-26) 

<„, = (4, • a* + 4, • a\ + At )• ^ • £,, , V* <= K, , Vy e N 0/ (4-27) 

K, + 4, + 4, = / , V* e K, , Vy e N D/ (4-28a) 

^,^,^£{0,7} , V*eK, , V/eN D/ (4-28b) 

x£ u ^,Sx- , V(/,y)eA (4-29) 

where 

Ar = subset of A for reservoir carry-over storage arcs; 

N/{/ = subset of N for reservoirs that have rule curves; 

No/ = subset of N for demand nodes subject to hedging rules; 

K, (or K^ ) = subset of reservoirs that supply demand nodey (or node //); 

D Jt = planned demand for demand node j; 

NK = number of reservoirs in K,; 
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PEj= percentage of the total water supply from reservoir k to demand node j, 

NK 

IP£>1; 

%j , £j , and X^j = binary (0-1) variables for reservoir k, 

S k mmi = minimum storage of reservoir k, 

S k ftrmj = firm storage of reservoir h, 

SLgaj = target storage of reservoir kr, 

S k max,, = maximum storage of reservoir t, 

a, and a\ = rationing factors for reservoir k, and all other variables are as defined before. 
Eqs. (4-22) ~ (4-28) ensure that the rule curves and hedging rules are observed in 
a multi-reservoir setting. Eq. (4-29) specifies the upper and lower bounds of each flow 
variable; for example this equation can represent the link capacity limitation, water 
supply limitation or reservoir storage limitation. The nonnegative decision variables 
include the flow and binary variables. 

It is important to note that the determination of k and K/ depends on the particular 
demand node j under consideration. For example, if demand node 1 can receive water 
from reservoirs 1 and 2, then k varies from 1 to 2, and K/ contains {R 1 , R 2 } for demand 
node 1 . In this example, we have PE\ + PE* = 1 . Likewise, if demand node 2 can receive 
water from reservoirs 3, 4 and 5; then, k =1 stands for reservoir 3, A: =2 stands for 
reservoir 4, and k =3 stands for reservoir 5, and K^ contains {R 3 , R 4 , R 5 } for demand 
node 2. In this example, it is obvious that PE\ + PE\ + PEl = 1 . 
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4.4 Solution Methodology 

There are several approaches for solving an MILP problem, but the most popular 
approach uses the branch-and-bound algorithm. With the branch-and-bound algorithm, 
the original problem is recursively decomposed into smaller problems to construct a tree 
of subproblems. These subproblems are systematically explored with the feasible integer 
solution and corresponding objective bounds to reach the optimal solution. The details of 
the branch-and-bound algorithm can be found in Nemhauser and Wolsey (1988) and 
Wolsey(1998). 

4.5 Numerical Example 

The following two-reservoir hypothetical simplified water distribution system is 
used to demonstrate the applicability and capability of the presented MILP model. 

System Description 

Figure 4-3 shows the configuration of this simplified water supply system. In this 
system, there are two supply sources, two reservoirs (RE1 and RE2), two demand nodes 
(DD1 and DD2), two diversion nodes, two terminal nodes and nine arcs. The water from 
the upstream sources flow into the reservoirs, and then the water will be stored or 
released to supply downstream demand nodes. Thus reservoirs play very important roles 
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in operating this system in order to determine the optimal water supply. To consider 
multi-time-period operation, each reservoir carry-over arc connects two adjacent time 
periods in a planning horizon to construct a large-scale water distribution network 
structure. 

Assuming that the planning horizon has 1 2 time periods. The inflow data of each 
inflow node used is shown in Table 4-1. The maximum storage of RE1 is 50 units and the 
maximum storage of RE2 is 30 units. The initial storage of RE1 is 35 units and the initial 
storage of RE2 is 25 units. Each reservoir has two corresponding rule curves and 
rationing factors (Table 4-2 and 4-3, and Figure 4-4). The rule curves are on a 12-time- 
period basis. The water demand of DD1 is 10 units and that of DD2 is 20 units in each 
time period. The maximum capacity of each arc is 100 units. Therefore, the capacity 
limitation is not a restricting factor for water delivery. 

The interpretation of the rule curves and of the hedging rules for reservoir RE1 is 
as follows: 

1 . If the beginning storage is greater than that of the target curve, the planned demand at 
each demand node is met at 100%. 

2. If the beginning storage is greater than that of the firm curve and less than that of the 
target curve, the planned demand for the demand node is met at 70%. 

3. If the beginning storage is less than that of the firm curve, the planned demand for the 
demand node is met at 50%. 
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The interpretation of the rule curves and of the hedging rules for reservoir RE2 is 
as follows: 

1 . If the beginning storage is greater than that of the target curve, the planned demand at 
each demand node is met at 100%. 

2. If the beginning storage is greater than that of the firm curve and less than that of the 
target curve, the planned demand at each demand node is met at 80%. 

3. If the beginning storage is less than that of the firm curve, the planned demand for 
each demand node is met at 60%. 

In the network configuration, demand node DD1 can receive water from reservoir 
RE1 and demand node DD2 can receive water from reservoir RE1 and RE2. Water 
supplies to these two demand nodes are required to satisfy the hedging rules in the 
optimization. Since DD2 can obtain waters from RE1 and RE2, in this case, assume 50% 
of total water supplying DD2 is required from RE1 and 50% is from RE2. Therefore, 
DD2 can be divided into two demand nodes, DD2 1 and DD2 2 with 10 units of water as 
planned demand respectively (Figure 4-5), then the above-mentioned overall MILP 
model can be implemented to optimize the system operation with hedging requirements. 

MILP Model 

The two hypothetical competing objectives considered in this case are: (1) 
maximization of water supply for the demand nodes; and (2) maximization of reservoir 
storage. The weighting method is used to combine the two objectives in the objective 
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function to reflect their relative importance. The weighting factors of water supply and 
reservoir storage are assumed 1 and 0.0001, respectively. In this situation, the priority of 
water supply is higher than that of the reservoir storage function. 

The maximization optimization problem can be converted to a standard 
minimum-cost problem by multiplying the objective function with a minus sign. 
Consequently, the 1 2-time-period MILP model for the simplified system can be 
formulated as follows: 



Mi "t " Z(C D ),,*,,,,- Z(cJ,,vi ( 4 " 3 °) 



i=l ieS DI jeN R , 

. (i.;)«A ('.y)eA R 

s.t. *(,,), =^„ ' V -/"eN, (4-31) 

Z *u,» ~ I x o,).. = , Vy e {N c uN fl J (4-32) 

0.')eA (/.;)eA 

I *(„>, " I x (,j), + J- = , Vy € N fl (4-33) 

(v.')eA (i.;)eA 

? ' = *(..,).,♦! -*(,.,)., < V ('".y) 6 A R ' V/' £ N K ( 4 ^4) 

£^. (4-22H4-29) 
where 

N/ = subset of N for inflow nodes; 
Ng = subset of N for junction nodes; 
Nd/v = subset of N for diversion nodes; 
INjj = amount of inflow to inflow nodey; 
Co = unit cost of supply to demand node; 
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Cr= unit cost of storage arc of the reservoir. 

The multi-objective function is formulated by Eq. (4-30). The continuity for the 
inflow, junction, diversion and reservoir nodes are denoted by Eqs. (4-31)~(4-34), 
respectively. In the optimization model, there are 344 decision variables, including 66 
integer variables and 370 constraints. Because the initial storage of RE1 and RE2 is 
assumed in zone 2, the integer variables for the beginning of the first time period does not 
need to be evaluated and the corresponding water supply for DD1 and DD2 is already 
pre-determined in the optimization model. 

In this study, the MILP optimization problem is solved with a commercial 
optimization program, LINGO. LINGO is capable of solving linear, nonlinear and integer 
programming problems. LINGO uses the branch-and-bound algorithm to deal with the 
integer variables. The details of LINGO are described in the LINGO 6.0 user's guide 
(1999). 

Using LINGO to solve this MILP model, the optimization results are shown in 
Table 4-4. The objective function value of the MILP model is -266.0. The results 
indicate that the reservoir storage and water supply simultaneously satisfy the rule curves 
and the hedging rules. For example, the solution shows that the beginning storage of 
reservoir RE1 at the second time period is in zone 2, and the corresponding optimal water 
supply for DD1 and DD2 1 is 7 units and thus meets the second hedging rule. 
Consequently, the capability of the proposed MILP model in handling hedging rules for 
reservoir operation in a water distribution system is verified by the above illustrative 
example. In addition, the results demonstrate that the main function of hedging for 
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reservoir operation is to utilize water more efficiently during a drought (i.e. conditionally 
restrain the current water uses for future water supply) and thus reduce the impact of 
drought to water users. 
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Table 4- 1 Water Sources (Unit) of the Simplified System 



Source 


Time 


Perioc 




1 


2 


3 


4 


5 


6 


7 


8 


9 


10 


11 


12 


1 

2 


10 
8 


10 
8 


10 
8 


12 
10 


12 
10 


12 
10 


14 
12 


14 
12 


14 
12 


8 
10 


8 
10 


8 
10 



Table 4-2 Rule Curves (Unit) Used in the Simplified System 



Time 


Rl 


El 


RE2 


Firm Curve 


Target Curve 


Firm Curve 


Target Curve 


1 


30 


50 


20 


30 


2 


27 


47 


17 


27 


3 


24 


44 


14 


24 


4 


21 


41 


11 


21 


5 


18 


38 


8 


18 


6 


15 


35 


5 


15 


7 


15 


35 


5 


15 


8 


18 


38 


8 


18 


9 


21 


41 


11 


21 


10 


24 


44 


14 


24 


11 


27 


47 


17 


27 


12 


30 


50 


20 


30 



Table 4-3 Hedging Rules Used in the Simplified System 



Reservoir Storage 



Rationin 



£ 



RE1 



Factor 



RE2 



Minimum Storage < S t < Firm Storage 

Firm Storage < S t < Target Storage 
Target Storage < S t < Maximum Storage 



50% 
70% 
100% 



60% 
80% 
100% 



S t : Beginning storage of RE1 (or RE2) in time period 



t. 
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Table 4-4 The MILP Model Results of the Simplified System 



Time 


Water Supply 


Reservoir Beginning Storage 


DD1 


DD2 


RE1 


RE2 


DD2' 


DD2^ 


Storage 


Zone 


Storage 


Zone 


1 


7 


7 


8 


35 


2 


25 


2 


2 


7 


7 


8 


31 


2 


25 


2 


3 


7 


7 


8 


27 


2 


23 


2 


4 


7 


7 


10 


23 


2 


24 


3 


5 


7 


7 


10 


21 


2 


24 


3 


6 


7 


7 


10 


19 


2 


24 


3 


7 


7 


7 


10 


17 


2 


24 


3 


8 


5 


5 


10 


17 


1 


26 


3 


9 


5 


5 


10 


20 


1 


28 


3 


10 


7 


7 


10 


25 


2 


30 


3 


11 


5 


5 


10 


19 


1 


30 


3 


12 


5 


5 


10 


17 


1 


30 


3 
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Figure 4-1 Rule Curves and Hedging Rules Represented by a Step Function 
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Figure 4-2 Rule Curves and Hedging Rules for a Multi-reservoir System 
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Figure 4-3 Schematic of Simplified System 
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Figure 4-4 The Rule Curves Used in the Simplified System 
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Figure 4-5 Configuration Transformation for the Simplified System 
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CHAPTER 5 MIXED-INTEGER NONLINEAR 
PROGRAMMING MODEL 



In general, the reservoir rule curves and hedging rules are pre-determined in the 
planning stage. However, during the reservoir operation, the sediment gradually changes 
the water storage of a reservoir and the population and economic growth in the reservoir 
serving area usually increases demands on the water supply. As a result, it is necessary to 
evaluate and modify the current reservoir rule curves and hedging rules to improve 
reservoir operation performance in a regional water distribution system. 

Chen (1995) developed a methodology incorporating genetic algorithms to 
modify rule curves for a single reservoir system. The results showed that the shape of rule 
curves was in oscillatory distribution and difficult to follow for real practice. Shih and 
ReVelle (1995) proposed a mixed integer-programming model to determine the optimal 
volumes of reservoir storage plus anticipated inflow that triggered several phases of 
hedging for reservoir operation during drought periods. Their model was applied for a 
single reservoir system. Neelakantan and Pundarikanthan (2000) presented a neural 
network-based simulation-optimization methodology to find the optimal reservoir storage 
levels for hedging. Although the case study of their research comprised of several 
reservoirs, the storages of these reservoirs were lumped together in the mathematical 
model. All the studies above have only considered obtaining the optimal reservoir 
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storages in order to use hedging rules for single reservoir systems. The rationing factors 
were pre-determined as given values in these models. 

This present study focuses on modifying both the current rules curves and 
hedging rules for multi-reservoir systems. Preferably, actual reservoir operation will 
easily follow the shape of the new rule curves which are obtained. 

Examining the overall MILP model described in Chapter 4, Eqs. (4-24H4-27) 
are linear constraints with given rule curves and rationing factors, and the binary variable 
set is one subset of decision variables. If the reservoir rule curves and hedging rules need 
to be evaluated, the rule curves and rationing factors become unknown and these linear 
constraints change to bilinear constraints, which are nonlinear. The nonlinearity is caused 
by the multiplication of one rule curve variable (or rationing factor variable) and one 
binary variable. Thus the original MILP model will change to a mixed-integer nonlinear 
programming (MINLP) model, which is a nonconvex optimization problem and is 
difficult to resolve. 

The general formulation of a MINLP problem with binary variables can be posed 
as(Floudas, 1995): 

Min. f(x,y) (5-1) 

s.t. h{x,y) = (5-2) 

g{x,y)<0 (5-3) 

x e X c R" (5-4) 

^F = {0,lf (5-5) 
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where 

x = a vector of n continuous variables; 

y = a vector of q binary (0-1) variables; 

fixy) = objective function; 

h(xy) = equality constraints; 

g(xj>) = inequality constraints. 

A MINLP problem is a NP-complete problem (Nemhauser and Wolsey, 1988) 
which consists of the continuous domain and the combinatorial domain. The presence of 
these two domains causes the complexity of this problem and creates the difficulties in 
solving it. To determine the global optimum of the MINLP problem in the nonconvex 
domain, which may consists of a number of local optimums, is NP-hard (Murty and 
Kabadi, 1987). 

5.1 Solution Methodologies for Solving MINLP Models 

To handle the MINLP problems, several types of methods had been reported in 
the past. These methods are: (1) combination method; (2) decomposition methods; (3) 
approximation methods; (4) heuristic combinatorial optimization methods; and (5) other 
methods. The following will investigate the mathematical characteristics of these 
methods and illustrate a number of application examples. 
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Combination Method 

As to solving a MINLP model, the intuitive manner is to use a combination 
methodology, which consists of a nonlinear programming (NLP) method and the branch- 
and-bound method. The basic idea of this combination method is to solve one NLP 
problem in each branch of the "tree" and proceed with this process until the final 
optimum is reached. In general, it is difficult to obtain the optimum by implementing this 
methodology because numerous branches along with their corresponding NLP problems 
need to be evaluated to finally achieve the optimal solution. If the dimension of one NLP 
problem or the integer variable set is large, the efficiency of this methodology is doubtful. 
In addition, the optimal solution obtained from this methodology can not be guaranteed 
as a global optimum at the nonconvex domain. 

Kim and Mays (1994) proposed an MINLP model for minimizing the cost to 
rehabilitate and replace pipes in a water distribution system. This model was solved by a 
combination scheme involving a branch-and-bound algorithm and a generalized reduced 
gradient method. In their research, for cases with large number (say >1000) of possible 
combinations, 1000 possibilities were randomly selected to represent the whole search 
domain. Because of the nonlinear characteristic and the possibility of numerous branch 
searching, this methodology may not be suitable to solve large-scale MINLP problems 
and the solution may be the local optimum. 
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Decomposition Methods 

The fundamental idea of the decomposition methods is to seek the upper bounds 
and lower bounds of the objective function value of a MINLP optimization problem. The 
optimal solution can be obtained by the convergence of upper bounds and lower bounds. 
The original problem is divided into two subproblems: the primal problem and the master 
problem. The solution information of the primal problem can provide the nonincreasing 
upper bounds. The solution information of the master problem can contribute to finding 
the nondecreasing lower bounds. Theoretically, with a limited number of iterations 
between solving the primal problems and master problems, the global optimal solution 
can be achieved by the convergence of the upper bounds and lower bounds. 

Two popular decomposition methods used to solve MINLP optimization 
problems are the Generalized Benders Decomposition (GBD) method and the Outer 
Approximation (OA) method. Floudas (1995) reported an excellent introduction 
describing GBD (and variants) and OA (and variants) for dealing with MINLP problems. 
The introductions of GBD and OA will be briefly described in the following sections. 

Generalize d Benders Decomposition (GBD) Method 

Benders (1962) developed a decomposition method for the solution of mixed- 
variables programming problems. In 1972, Geoffrion generalized Benders' 
decomposition approach to obtain the optimal solution for certain nonconvex NLP and 
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MINLP problems. In the GBD algorithm, the upper bound of optimal solution and 
Lagrange multipliers associated with equality and inequality constraints are obtained 
from the primal problem in which the complicating variables are fixed. By fixing the 
complicating variables, it becomes easier to solve the original problem. For instance, the 
integer variables in a MINLP problem can be the complicating variables. The master 
problem is developed based on the nonlinear duality theory involving the Lagrange 
multipliers that are resulted from the primal problem. The new set of complicating 
variables is obtained after solving the master problem and used in the next primal 
problem as fixed variables. Proceeding with the iteration between the primal problem and 
master problem the optimal solution can be achieved by the convergence of upper bounds 
and lower bounds. Bagajewicz and Manousiouthakis (1991) reported an analysis on 
various implementations of GBD and arrived at several conclusions. 

Jacobs et al. (1995) reported a stochastic optimization model to deal with the 
optimal scheduling of hydropower generation problem for a multi-reservoir system. The 
model was conducted based on network formulation and solved by GBD. The authors 
found that the performance of GBD in computing efficiency was very encouraging in 
solving this kind of problem. 

Cai et al. (2001) applied GBD to solving two large nonconvex water resources 
optimization problems. In their model formulation, the elastic slack variables were 
introduced and incorporated in the objective function as penalty terms to guarantee the 
feasibility of the primal problem. Compared with two NLP solvers, MINOS (Murtagh 
and Saunders, 1982) and CONOPT (Drud, 1994), the authors found that GBD could 
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obtain similar objective function value in less run time if the complicating variables were 
carefully selected. 

Outer Approximation ( OA) Method 

The concept of OA was originally proposed by Gomory(1958, 1960), Cheney and 
Goldstein(1959), and Kelley (1960). Duran and Grossmann (1986a,b) presented the OA 
algorithm to solve certain MINLP problems. The basic procedure of OA is similar to 
GBD except for the derivation of the master problem. In the OA algorithm, the master 
problem is derived using the solution of the primal problem and outer linearizations of 
the nonlinear functions around the solution of the primal problem. 

Karatzas and Pinder (1993) presented a methodology in which the OA method 
was incorporated to deal with the groundwater quantity management problem. Since the 
fixed charges were formulated as an exponential term in the objective function, the 
resulting problem was a concave minimization problem. Two examples were presented to 
demonstrate the efficiency and robustness of this methodology in finding the global 
optimal solution. 

Watkins and Mc Kinney (1998) used GBD and OA algorithms to solve special 
water resources MINLP problems in which discrete and nonlinear terms were in the cost 
functions, and then compared the performance of these two algorithms. Although it was 
difficult to claim which one of these two algorithms was better, the authors suggested that 
GBD may be better suited for dealing with larger scale planning problems. 
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Approximation Methods 

The primary difficulty for solving a MINLP problem is due to the combinatorial 
characteristic caused by the discontinuous/integer variables. One way to eliminate this 
obstacle is to use continuous variables to approximate the discontinuous variables to 
obtain a continuous search domain. Thus the original MINLP model can be converted 
into a nonlinear programming (NLP) problem, which can be solved by conventional 
nonlinear programming algorithms. 

Two approximation methods proposed in the water resources field, especially in 
the groundwater management area, are: (1) the penalty coefficient (PC) method, which 
includes the polynomial penalty coefficient (PPC) method and the exponential penalty 
coefficient (EPC) method, and (2) the pseudo-integer (PI) method. Note that using these 
approximation methods may create a complex nonlinearity in the resulting NLP problem 
and the optimal solution may be the local optimum, but is not assured to be the global 
optimum. 

Mckinney and Lin (1995) reported a groundwater MINLP model to deal with the 
optimal aquifer remediation design problem. In their formulation, the binary variables 
were involved in the objective function. They used PPC method to approximate the 
binary variables to convert the original MINLP problem into a NLP problem. Compared 
with PI and EPC method, the authors claimed that the model performance of PPC was 
better than PI and EPC in solving MINLP optimal remediation design problems. 
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Heuristic Combinatorial Optimization Methods 

In an optimization problem, if the discrete decision variables are involved and a 
solution is searched in a finite set of possible solutions to meet the optimum, this 
optimization problem is a combinatorial optimization problem. Searching for the global 
optimal solution of a combinatorial optimization problem is a challenging task. As a 
result, a number of heuristic search methods, such as Genetic Algorithms (GA) and 
Simulated Annealing (SA), are proposed to solve combinatorial optimization problems. 
What these methods have in common is the lack of a gradient mathematical 
characteristic. The following will briefly summarize the concept of GA and SA and 
report some example applications. 

Genetic Al gorithms (GA } 

GA was originally developed by Holland in 1975. Motivated from biological 
evolution, GA uses a set of analogous process and operators: selection, crossover, 
mutation, reproduction, and replacement, to search for a global or near-global optimal 
solution of a complex optimization problem. Unlike the traditional gradient-based 
algorithms, the derivatives are not required in GA. A population of starting chromosomes 
that are randomly generated to represent a set of initial solutions goes through the GA 
process and operations and improves the objective value generation by generation. By 
evolving chromosomes without the necessity of gradient evaluations, GA has the ability 
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to jump out of local optima and achieve the global optimum. Readers can see Goldberg 
(1989) for a comprehensive introduction of GA. 

Ritzel and Eheart (1994) presented the application of GAs on a multiple objective 
groundwater pollution containment problem. Finding the optimal solutions on the trade- 
off curve between the reliability and cost of a hydraulic containment system was the main 
task in their research. Compared with the mixed-integer chance constrained 
programming, the results obtained from their example study showed that GA could 
obtain similar or better optimal solution with less computer time. 

McKinney and Lin (1994) reported a methodology in which groundwater models 
and GAs were combined to solve three complicated groundwater management problems: 
maximum pumping from an aquifer, minimum cost of water supply development, and 
minimum cost of aquifer remediation. The application results showed that the 
performance of GA was as good as or superior to that obtained by LP and NLP. In 
addition, GA was verified as having the ability to efficiently search for the global or near- 
global optimal solution of those groundwater management and design problems. 

Oliveira and Loucks (1997) presented a methodology incorporating GA to derive 
multireservoir operating policies. The operating policies can indicate the total release of 
reservoirs, the individual reservoir storage targets, and the existing total storage volume 
in the system. In their approach, a simulation model was used to evaluate each generated 
operating policy in GA. The process of updated policy generation and evaluation 
proceeded until the stopping criterion was met. The applicability of this methodology was 
demonstrated by two hypothetical reservoir systems. 
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For more on water resources GA applications, readers are referred to Goldberg 
and Kuo (1987), Simpson et al. (1994), Halhal et al. (1997), Savic and Walters (1997), 
Montesinos et al. (1999), and Cai et al. (2001). 

Simulated Annealinv (SA) 

SA was first proposed by Kirkpatrick et al. (1983). This algorithm was inspired 
from a physical process of thermodynamics, in which a substance is heated and then 
slowly cooled to obtain a strong crystalline structure. In SA, a new state of a system is 
constructed by random displacement. During the cooling procedure, if the energy of the 
current state is lower than that of the previous state, the change is accepted 
unconditionally and the system is updated. If the energy of the current state is greater 
than that of the previous state, the new state is accepted probabilistically by the 
Metropolis criterion. Thus SA has the ability to reach the global optimum or near-global 
optimum instead of trapping in the local optimum. 

Dougherty and Marryott (1991) reported the application of SA on the 
optimization of groundwater management problems. The mathematical characteristics of 
SA were examined and investigated. By four illustrative applications of SA, the authors 
claimed that SA was flexible and had the potential for solving combinatorial groundwater 
management problems. 

Unlike the conventional optimization methods used in the past (e.g., linear 
programming, nonlinear programming, and dynamic programming), a SA-based 
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approach was proposed by Cunha and Sousa (1999) to deal with the optimal water 
distribution design problem. Compared with the solutions of past well-known studies, the 
authors claimed that their approach could provide solutions with high quality for network 
design problems. 

For more on water resources SA applications, readers are referred to Marryott et 
al. (1993), Mauldon et al. (1993), Loganathan et al. (1995), Zheng and Wang (1996), and 
Skaggsetal. (2001). 

Other Methods 

The other methods include several heuristic approaches for some specific 
problems. In literature, it is very difficult to find a robust methodology, which is 
universally accepted to solve MINLP problems. If there is no existing methodology that 
can efficiently solve a MINLP problem, a heuristic methodology is usually developed 
based on the special mathematical structure of this problem. Therefore, in general, one 
heuristic approach is only suitable to one kind of problem; this approach may not be used 
to deal with other problems with different mathematical structures. 

Hsu and Yeh (1989) reported a heuristic approach to solve a MINLP problem for 
experimental design problems in groundwater hydrology. This approach was an implicit 
trial-and-error method based on experience to pre-choose the pumping wells before the 
optimization calculation. This approach took advantage of the mathematical 
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characteristics of their model formulation. Although this approach can efficiently solve 
their MINLP problem, the optimal solution was not assured as the global optimum. 

Kwanyuen and Fontane (1998) presented a groundwater planning problem which 
considered the minimization of the installation and operation cost and was formulated as 
a MINLP problem using the response matrix method. In their research, a heuristic 
branch-and-bound method was developed based on taking advantage of the structure of 
the groundwater optimization problem by limiting the search process to the high potential 
branching only. Compared with other methods, the authors claimed that this methodology 
was more robust in finding the global optimal solution. 

In this study, GBD, PC, PI, and SA are chosen to solve the presented MINLP 
model to achieve a new set of reservoir rule curves and hedging rules. In addition, a 
transformation technique will be proposed to solve this MINLP model as well. The 
results of these methods will be investigated and compared. The description of the 
application of these methods will be described in the following sections. 

5.2 Generalized Benders Decomposition 

In the GBD algorithm, a MINLP problem is divided into two subproblems: the 
primal problem and master problem. The idea of the primal problem is to fix the 
complicating variables to change the original MINLP problem into an easier problem. 
The primal problem can be formulated as: 
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Min. f(x,y m ) (5-6) 

x 

s.t. h(x,y m ) = (5-7) 

g(x,y m )<0 (5-8) 

xeXcR" (5-9) 

In the primal problem, the set of complicating variables, y, is fixed in ot* 

iteration, /". Thus x is the only decision variable set in the primal problem and can be 

solved without much effort. However, since not only the feasible solution, but also the 

infeasible solution may occur in the primal problem, this kind of consequence will affect 

the derivation of the master problem. 

The derivation of the master problem is based on the nonlinear duality theory. 
After a sequence of mathematical derivation, the master problem is transformed into a 
relaxed master problem. The detailed discussion of the master problem is ignored here. 
Readers are referred to Floudas (1995) for details. The relaxed master problem in GBD is 
formulated as: 

Min. Mb (5-10) 

s.t. M B >£(yU k ,M"),k = \,2 K (5-11) 

>£(>>; A V), / = 1,2,...^I (5-12) 

where 
H B = a scalar; 

k = optimal multiplier vectors for the equality constraints at A** 1 feasible primal problem; 
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//*= optimal multiplier vectors for the inequality constraints at k^ feasible primal 
problem; 

k = optimal multiplier vectors for the equality constraints at f* 1 infeasible primal 
problem; 

H = optimal multiplier vectors for the inequality constraints at / ch infeasible primal 

problem; 
£(y; A*,//* ) = support function, derived from the solution ofk A feasible primal problem; 

l(y;^ <M) = support function, derived from the solution of 7 th infeasible primal problem. 

As the procedure proceeds, the size of the relaxed master problem becomes larger 
and larger because more and more constraints associated with the solution information 
from the previous feasible or infeasible primal problems are added in this problem. If the 
optimal solution is not achieved in a small number of iterations, the efficiency of solving 
the resulting relaxed master problems will be doubtful. 

To solve the above-mentioned MINLP problem to update the current rule curves 
and hedging rules, the variant 3 of GBD, v3-GBD (Floudas, 1995), which is called global 
optimum search (GOS) algorithm (Floudas et al., 1989), is implemented. 

Different GBD variants have corresponding relaxed master problem formulation. 
For v3-GBD, the relaxed master problem is: 

Min. n B (5-13) 

y*Y fa 

s.t. n B >L{x\y,X k ,n k ),k = \,2^..,K (5-14) 
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0>L(x l ,y,A.' ,n'), / = U n ..^ (5-15) 



where 



L(x k ,y,A k , M k ) = f(x k ,y) + Z kT h(x k ,y) + M kT g(x k ,y) 
ZCx.yj'.^V/^x'.yJ + ^'^x'.y) 

are the Lagrange functions associated with the optimal solution (x* or x 1 ) of the primal 
problem. 

The flowchart of v3-GBD is shown in Figure 5-1. The detailed steps of v3-GBD 
are: 
Step 1: Select an initial^ 1 with which the first primal problem is feasible. Solve the first 

primal problem with y l and obtain an optimal solution x 1 and multiplier vectors 

X and u 1 . Set the counters k = 1, / = 0. Let the optimal objective function value 

beX^'iV 1 )- Setyfjc 1 ^ 1 ) as the current upper bound, UBD. Select a nonnegative 

convergence tolerance, e. 
Step 2: Solve the relaxed master problem to obtain the optimal solution y and// B . Set 

Mb as the current lower bound. LBD. If UBD - LBD < e, then terminate the 

procedure. 
Step 3: Solve the primal problem with the new y set, which is y . Then the optimal 

solution can be either feasible or infeasible: 
Step 3 a: If the primal problem is feasible, set k = k + / and store the values of multipliers, 

k and //. Let the optimal objective function value be fix*, y ). Update the upper 



52 



Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 



bound UBD = min{UBD,/x*, y')}. If UBD - LBD < e, then stop. Otherwise, 
return to step 2 with jc*, A* and /A 
Step 3b: If the primal problem is infeasible, set / = / + 1. Solve the following relaxed 
primal problem (Floudas et ah, 1989) to obtain an optimal solution x 1 and the 

multipliers k and /i . Then return to step 2. The formulation of the relaxed 



primal problem is: 



Min. a (5-16) 

xeA' 

s.t. h(x,y')-a<0 (5-17) 

-h(x,y')-a<0 (5-18) 

g(x,y)-a<0 (5-19) 

a > (5-20) 



where 



a = nonnegative slack variable. If a is zero, the infeasibility of the primal problem is 
withdrawn. 

The implementation of v3-GBD will be applied to a hypothetical simplified 
system and the discussion will be described later. 

S3 Penalty Coefficient Methods and Pseudo-integer Method 

The difficulty of solving the proposed MINLP model is due to the existence of 
binary variables, which are discontinuous variables. As a result, the main purpose of 
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implementing approximation methods is to use continuous variables to approximate 
binary variables. Thus, the original MINLP problem is converted into a NLP problem, 
which can be solved by conventional nonlinear programming algorithms. 

Polynomial Penalty Coeff icient (PPC) Method 

The basic idea of PPC (McKinney and Lin, 1995) is to use a continuous variable 
and a division form to approximate the original binary variable. Let 7* , represent X k wt . 
Implementing PPC to approximate the binary variables of the presented MINLP model, 
7*, can be formulated using x" y) , : 

Vwj =— (5-21) 

x u.ju + P 

where w= 1 ,2, and 3 and p is a very small positive value, 0<p« 1 . Here, x" U)J represents 

the reservoir storage variable. Introducing this formulation into Eqs. (4-24)~(4-28a), the 
original MINLP model will be converted to a continuous NLP model with a high degree 
of nonlinearity, which is caused by the multiplication terms in which the division forms 
are involved. 

The original mathematical characteristics of the three binary variables, which 
specify that only one of the three zones in the reservoir can be active, can still exist by 
PPC. For example, let w = 1, if x\ lJU = 0, then 7*,= and A{, = 0; if x\ i])t * 0, 7*, 
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and X\ x can be approximately considered as 1, and ^ = 0, X\ x = 0, and xf U)J = 0, x (iJU 

= 0. 

The reason of using the reservoir storage variables to approximate the integer 
variables is to represent the presence of which one of the three reservoir storage variables 
and the corresponding binary variables should be active. 

Exponential Penalty Co efficient (EPC) Method 

The concept of EPC (Karatzas and Pinder, 1993) is similar to PPC. A continuous 
variable and an exponential term are used to approximate the original binary variable. In 
this study, when EPC is implemented to approximate the binary variables, let 7* , 
represent X k WJ , and 7*, can be formulated using x" J)t : 

7 * r= y- e " x <'^ (5-22) 

where w = 1,2, and 3 and u is a very large positive value. Here, x ( 7 y)/ represents the 

reservoir storage variable. As a result, the original MINLP model can be transformed to a 
continuous NLP model with a high degree of nonlinearity, which is mainly caused by the 
multiplication terms in which the exponential forms are involved. 

The reason for using the reservoir storage variable in EPC is the same as PPC. 
The original mathematical characteristics of the three binary variables, which represent 
that only one of the three zones in the reservoir can be active, can still exist by EPC. For 
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example, let w = 1, if x x UJ)J = 0, then 77*,= and 4*,= 0; if x\ lj)t * 0, 77*, and Af, can 
be approximately considered as 1, and ^ = 0, A% a = 0, and xf IJU = 0, x ( 3 /y) , = 0. 

Pseudo-int eeer (PI) Method 

The concept of PI (Lall and Santini, 1989) is to introduce a continuous variable 
and an inequality constraint to replace the original binary variable. For the original 
MINLP model in this study, let 77*, represent A*, and introduce the following 
constraint: 

nirU-ni,)*** ( 5 - 23 ) 

where 77*, is a continuous variable, < rj k WJ < 1 . If rj k WJ = 0, means A.^ = in the original 
model; if 7*,= 1. represents A* x = 1. Therefore, the original MINLP problem becomes a 

NLP problem with continuos variables and constraints. Note that it may be difficult to 
deal with Eq. (5-23) by the conventional gradient-based method when it is preferred to 
obtain the 0-1 solution. 

The main difference between PC and PI methods is that in the PC methods the 
k WJ variables are computed as functions of the reservoir storage variables, while in the PI 
method the A k WJ variables are computed by the introducing continuous variables and 
constraints. The results of PC are either zeros or values which are almost equal to one, 
and the results of PI are enforced to be zeros or ones. 
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5.4 Simulated Annealing 

Annealing is a physical process in which a solid is heated up and then slowly 
cooled down to reach a state with the lowest energy. At the initial high temperature, the 
molecules in the material have more freedom with high energy to reach different states. 
As the temperature is declined, the energy is decreased as well. If the material is cooled 
properly, the molecules will reach an ordered crystalline configuration with minimum 
energy. If the material is cooled rapidly, there will be defects and glass-like intrusions 
created inside the material. Thus, this material does not reach the state with minimum 
energy and ends in a polycrystalline state with higher energy. The physical scheme of SA 
is shown in Figure 5-2. 

The SA optimization method is considered analogous to annealing: the energy 
represents the objective function, the configurations of molecules represents the 
configurations of control parameter values and decision variables, and seeking the state 
with lowest energy represents finding the global optimal solution. 

Let s represent the solution set of a minimum optimization problem and c 
represent the objective function of 5. In SA, a new solution (5,+/) is randomly generated in 
the neighborhood of the current solution (s,). As to the random generation of a new 
solution, readers can see Corana et al. (1987) for details. The difference between s, and 
5,*i is that only one element of s, is modified. Deciding to accept s^i or not is judged by 
the Metropolis criterion: 

If ( c^i - c, ) < 0, then accept s i+ i; 
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• C i C i+1' 



else accept or reject s l+/ with the acceptance probability/?: p = exp(— — ) . 

Tk is the current temperature after k times of temperature reduction. A 

pseudorandom number, pr, is generated in the range of [0,1] and compared with 

p. If pr < p, then accept *,,./, otherwise reject it. 
If s,~i is accepted, then update the current solution with it; otherwise, keep the current 
solution. 

Examining the Metropolis criterion, it is apparent that 7* is very important in the 
acceptance judgment. A higher Tk causes a higher acceptance possibility for the new 
solution and it implies that both uphill and downhill moving in the search domain are 
possible. Conversely, a lower Tk creates a lower acceptance possibility for the new 
solution and it implies that only downhill search is preferred to appear in the 
minimization problem. Based on this criterion, SA has the ability to avoid the local 
optimum and reach the global or near-global optimum after numerous evaluations. 

A set of control parameters employed in SA is called an annealing schedule. 
Generally, the determination of these control parameters relies on a trial-and-error 
experiment to obtain an encouraging optimal solution. Among these parameters, the 
influential ones are: the initial temperature (T ), the maximum number of iteration under 
a temperature (Nr), and the temperature reduction rate (Rr)- In SA, at each temperature 
the simulation must proceed long enough to reach a steady state. This is called 
thermalization (thermal equilibrium). The parameter, Nr, which should be appropriately 
large, is used to reach thermalization. To represents a sufficient high temperature in the 
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beginning of the annealing process and Rr determines how slowly this system should be 
cooled. Due to these characteristics, these parameters play crucial rules on deciding the 
performance of model convergence and the computer running time in the SA process. 
The stopping criteria used in SA can be decided by the user, e.g. specify a maximum 
number of evaluations, or assign a numerical error tolerance for the optimal solution. 

A promising control parameter set in SA can guarantee a global optimization 
solution. However, developing an excellent parameter setting standard is a very 
challenging work. In literature, a number of studies are already presented to suggest the 
determination of SA parameters (Kirkpatrick, 1984; Aarts and Van Laarhoven, 1985; 
Huang et al., 1986; Corana et al., 1987; Johnson et al., 1989). Nevertheless, none of the 
conclusions of the studies above is guaranteed to be applicable to all SA implements. For 
a SA application, a systematic analysis to determine the efficient SA parameters is critical 
before a global optimum is claimed. 

The main difficulty of solving the proposed MINLP model is because the 
multiplication terms of the continuous and binary variables create a high level of 
nonlinearity. A hybrid SA approach is developed to implement SA in solving this model 
to obtain a new set of rule curves and hedging rules. 

The basic idea of this hybrid SA is to use SA to evaluate the complicating 
variables, which can be binary variables or rule curve variables, and then use an 
optimization solver to solve this model with given complicating variable solution. By 
continuing to proceed with the iteration between SA and solving this model with given 
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complicating solutions until the stopping criterion is met, the optimal solution can be 
achieved. 

In this study, LINGO is used as the optimization solver. Figure 5-3 shows the 
flowchart of the hybrid SA approach. The steps of the hybrid SA are: 
Step 1 : Select a control parameter set, which includes To, Nt, Rt, and initial solution so- 

Calculate the objective function of s by LINGO. Let / = 1 and k = 0. 
Step 2: Randomly generate a new solution, s„ and calculate the corresponding objective 

function by LINGO. 
Step 3: Evaluate the new solution by the Metropolis criterion to decide to accept this 

solution or not. And update the current solution. 
Step 4: Check if the total evaluation number under a constant temperature so far is equal 

to Nt or not. If so, reduce the temperature, the new temperature T k +, = R T -T K 

and let k = k+1, then go to step 5. If not, let / = /+/ and go to step 2 and resume 

this process. 
Step 5: Check if the stopping criteria are met or not. If so, terminate this process and 

consider that the best solution so far is the optimal solution. If not, go to step 2 

and continue this process. 

In this study, a SA program developed by Goffe (1994) is modified and linked 
with LINGO to develop the hybrid SA approach. This program is constructed based on 
the SA algorithm proposed by Corana et al. (1987). 
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5.5 Transformation 

In the original MINLP problem, a product term in which a continuous variable 
and a binary variable are involved causes the difficulty in solving this problem. To 
overcome this difficulty, unlike the above-mentioned MINLP solving methods, a 
transformation technique (Williams, 1999) is used herein to transform the proposed 
MINLP problem into a MILP problem to eliminate the nonlinearity. 

Assume x is a continuous variable and X is a binary variable. Let z = xk and it is 
apparent that if k = 0, z = 0, and conversely if A. = 1, z = x. The transformation 
formulations which are based on Williams (1999) are modified and stated as follows: 

2 = x mm -^ + y (5-24) 

y-{^-x mm )-^<0 (5-25) 



-x + y<-x mw 

x -y + ( x max -x mn y^<x n 

y>0 



(5-26) 
(5-27) 
(5-28) 



where 

Xmm = lower bound of x; 

Xmax = upper bound of x. 

If A. = 0, then z = which is the same as the original result. Conversely, if X. = 1 
then z = x. Although this transformation technique creates two more continuous variables 
and four more linear constraints, the original MINLP problem can be transformed into a 
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new MILP which can be solved by any standard integer programming algorithm, e.g. the 
branch-and-bound method. This is the main advantage of the transformation technique. 

5.6 Numerical Example 

The simplified system presented in Chapter 4 is used herein to achieve a new set 
of rule curves and hedging rules by solving the MINLP model via the solution 
methodologies presented above. The MINLP model for this case study consists of Eqs. 
(4-22) ~ (4-29) with unknown rule curve and rationing factor variables and Eqs. (4-30) ~ 
(4-34). The original MINLP problem has 408 variables (including 66 binary variables) 
and 372 constraints (including 120 linear constraints and 252 nonlinear constraints). 

Inherently, the domain to search for a new rule curve should start from the 
minimum reservoir storage to the maximum reservoir storage. However, according to a 
number of model tests, even though an encouraging objective function value can be 
achieved, this kind of search domain results in a set of new rule curves which are in 
randomly, oscillatory distribution. Thus the reservoir operators can not operate the 
reservoir based on the new rule curves because the shape of rule curves is difficult to 
follow. As a result, in this study, a set of lower bounds and upper bounds based on the 
current setting are assumed to obtain the new rule curves and hedging rules. Similarly, 
the range assumption is also applied to search for new rationing factors to avoid a large 
difference between cti and 0:2. 
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For this simplified system, the lower bound of each value of each rule curve is 
assumed as the original value minus 2 units. The upper bound of each value of each rule 
curve is assumed as the original value plus 2 units. If the upper bound is equal to or 
greater than the maximum storage of reservoir, the upper bound is set as the maximum 
storage. The lower bound of aj is assumed as the original value minus 0. 1 . The upper 
bound of 0C2 is assumed as the original value plus 0. 1 . The searching ranges for new rule 
curves and rationing factors are shown in Table 5-1 and 5-2. 

The implementations of the MINLP solution methodologies are described and the 
results are discussed in the following sections. In these methodologies, LINGO is used as 
a solver to handle an optimization problem. 

Generalized Benders Decomposition 

Solving the MINLP model by GBD, the binary variables are chosen as the 
complicating variables. As a result, the primal problem is an LP problem and the relaxed 
master problem is a MILP problem. In each primal problem, the continuous decision 
variables include the water flows in this network system, the reservoir rule curves and the 
rationing factors. In the relaxed master problem, only one continuous variable is involved 
and the other decision variables are binary variables. 

According to the GBD algorithm, the initial binary variable guess must allow the 
first primal problem to be feasible. After a number of tests, a good initial integer guess is 
achieved: for each reservoir, the storage at the first 7 time periods is specified in zone 2, 
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i.e. /lf,= 0, J%j= 1 and A^ a = 0, / = 1~7, and the storage at the rest of the time periods is 
specified in zone 1, i.e. ^= 1, >#,= and A% t = 0, r = 8-12. Here, the nonnegative 

convergence tolerance, e, is assumed as 0.01. 

In each primal problem, there are 336 continuous variables and 360 linear 

constraints. In each relaxed master problem, there is one continuous variable and 66 

binary variables. The number of constraints is varied because one more constraint is 

added after one primal problem is solved. 

After a number of iterations between the primal and relaxed master problems, this 

study found that the convergence of the upper bounds and lower bounds is extremely 

slow. This is mainly because: 

1 . From the mathematical viewpoint, in each iteration, the model would like to minimize 
the objective function value of the relaxed master problem to be as small as possible. 
Therefore, based on the formulation of the relaxed master problem, the optimization 
model would provide as many A£,= 1 solutions as possible. This situation implied 

that the model would prefer the reservoir storage solution of each time period to be in 
zone 3. If the sources were plentiful, this situation may not generate an infeasible 
primal problem. However, this situation created infeasibility in the primal problem 
when the sources were not plentiful. Because the main purpose of incorporating 
reservoir rule curves in reservoir operation is to provide operating guidelines for 
drought periods; thus, the A£,= 1 solution was not the preferable result. This is a 
conflict situation between the physical system and the mathematical model. 
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2. The upper bound values kept as constant, which was the objective function value of 
the first feasible primal problem. Examining the primal problem of each iteration, 
most of the primal problems were infeasible, i.e., most of the binary solution of the 
previous relaxed master problem caused the infeasibility in the next primal problem. 
Thus a number of relaxed primal problems with artificial slack variables needed to be 
solved. 

3. The order of the range of objective function values of both primal and relaxed master 
problem was 10 2 . Although the lower bound values increased in each iteration, the 
general improvement performance of lower bounds was in an order of 10" . In 
addition, a number of results showed that the lower bounds kept nonincreasing. This 
condition implied that the convergence of the upper bounds and lower bounds was 
extremely slow and that the number of iteration runs might be multitudinous. 

4. As the iteration procedure proceeded, the size of the relaxed mastered problem 
became larger and larger because one additional constraint, in which most of the 
binary variables were involved, was added after one primal problem was solved. Thus 
the solving time became longer as well to solve a relaxed master problem in each 
iteration. This kind of condition resulted in space and memory problems when 
running the model on a PC. 

Although GBD was demonstrated as a suitable methodology in dealing with 
nonconvex water resources management problems in literature (Watkins and McKinney, 
1998; Cai et al., 2001), in this study, unfortunately, it may not be applicable to solving 
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the presented MINLP model because of the poor convergence performance. Even though 
various initial binary solutions were tested, the convergence still did not perform well. 

Penalty Coefficient Methods and Pseudo-inteeer Method 

To solve the resulting NLP problems by EPC, PPC, and PI methods, LINGO uses 
the generalized reduced gradient (GRG) algorithm. GRG is a gradient-based algorithm 
that is frequently used to solve large-scale nonlinear programming problems. Details of 
the GRG method can be seen in Lasdon and Waren (1978), and Drud (1994). 

In the NLP model from EPC method, there are 336 continuous variables, 108 
linear constraints, and 276 nonlinear constraints. The number of nonlinear constraints is 
over 70% of total constraints. The u value used herein is 10 6 . 

The number of variables and constraints of PPC method is the same with EPC 
method. The p value used herein is 1 0" 6 . 

In the NLP model from PI method, there are 408 continuous variables, 216 linear 
constraints, and 240 nonlinear constraints. The number of nonlinear constraints is over 
50% of total constraints. 

Solving the resulting NLP problems, the objective function value of EPC, PPC, 
and PI are -219.9, -220.0, and -264.0, respectively. The results show that EPC and PPC 
have similar objective function value. In addition, PI performs better than EPC and PPC. 
Using these three methods to transform the original MINLP problem into three 
continuous NLP problems can eliminate the searching procedure in the discontinuous 
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domain. However, in the nonlinear constraints of the resulting NLP problems, the 
multiplication terms in which the division or exponential forms are involved cause a high 
level of nonlinearity. In other words, it is possible that numerous local optimum exist in 
the search domain. Consequently, it is very difficult to obtain the global optimum for the 
resulting NLP problems by traditional gradient-based algorithms. 

Since the task is to obtain new rule curves and hedging rules, the objective 
function values of these three models are preferred over the objective value of the current 
rule curves and hedging rules, which is -266.0. Unfortunately, the results of these 
methods did not improve the objective function. The reasons may include that too many 
local optimums exist or that the initial guess is not well selected. Thus the solution is 
easily trapped in a local optimum. 

Simulated Annealine 

Two kinds of hybrid SAs are proposed to solve the original MINLP model for this 
simplified system and to provide a comparison of the results. 

The first approach, HSA1, the binary variables (discontinuous variables) are 
chosen as the complicating variables evaluated by SA. With given binary variable 
solution, the resulting model becomes an LP problem. Thus solving the original MINLP 
problem becomes solving a number of linear programming problems with their 
corresponding binary variable solutions. In each LP problem of HSA1, there are 336 
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variables and 360 constraints. The SA control parameters used in this case are: To = 80, 
Ar r =20,and/f r = 0.5. 

The second approach, HSA2, the rule curves and rationing factors (continuous 
variables) are chose as the complicating variables evaluated by SA. With the given rule 
curves and rationing factors solution, the resulting model becomes an MILP problem. 
Thus solving the original MINLP problem becomes a matter of solving a number of 
MILP problems with their corresponding rule curves and rationing factors solutions. 
Solving the MINLP model by HSA2, there are 48 rule curves values and 4 rationing 
factors which need to be evaluated. In each MILP problem of HSA2, there are 356 
continuous variables, 66 integer variables, and 384 constraints. The SA control 
parameters used in this case are: T = 80, N T = 20, and R T = 0.5. 

Figure 5-4 shows the relationship of objective function value and the number of 
evaluations in HSA1 and HSA2 approach. The relationship represents the best objective 
function value found so far after a number of iterations. A comparison of HSA1 and 
HSA2 in solving the MINLP problems for this simplified system reveals the following: 

1. Because HSA1 and HSA2 are iteration approaches, a myriad of LP problems are 
solved in HSA1 and a myriad of MILP problems are solved in HSA2 to obtain the 
optimal solution. 

2. The numbers of the decision variables and constraints of HSA2 are larger than HSA1 . 
In addition, the binary variables are involved in the resulting MILP model in HSA2. 
Therefore, the computer solving time of HSA2 in each iteration is longer than that of 
HSA1. 
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3 . In HS A 1 , infeasibility may occur. This kind of solution is omitted in S A evaluation. 

4. The objective function value of HSA1 is -281.0 and that of HSA2 is -279.1. HSA1 
and HSA2 can obtain similar objective function values. 

Transformation 

Using the transformation technique, there are 650 continuous variables, 66 binary 
variables, and 1032 linear constraints in the resulting MILP model. The MILP model is 
then solved by LINO to obtain a new set of rule curves and hedging rules. The objective 
function value obtained is -28 1 .0, which is superior to the objective function value under 
the current rule curves and hedging rules. Consequently, the results suggest that the 
operation for this simplified system can reach better performance if the new rule curves 
and hedging rules are implemented. 

Overall Discussion 

To reevaluate the rule curves and hedging rules of the hypothetical simplified 
system, the above-mentioned methodologies are applied to solve the MINLP model and 
the results are analyzed and discussed above. The characteristics of these methods are 
summarized in Table 5-3. 

In these methodologies, the MINLP model solved by transformation technique 
can obtain the most encouraging objective function value. In other words, compared with 
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other methodologies, the transformation technique is more suitable for solving the 
MINLP model to reevaluate the current rule curves and hedging rules. The optimal 
results by the transformation technique are shown in Table 5-4. The new reservoir rule 
curves and hedging rules are shown in Table 5-5 and Table 5-6. 

Examining the optimal rule curve and rationing factor results, it is obvious that 
most of the values reach either the lower bounds or the upper bounds. This situation 
implies that there is potential to obtain a better objective function value if the search 
range becomes larger. However, with a larger search range, the new rule curves are likely 
to have oscillation distributions and be difficult to follow for real operation. 

To discover how the objective function responds on the search range, the 
sensitivity analysis is implemented by adding four more search range sets and the results 
are compared with the original setting (Table 5-7). The results of these five scenarios are 
shown in Table 5-8 and Figures 5-5-5-9. Note that the result of ai for RE2 in each 
scenario can be ignored. The reason is because the model results show that the reservoir 
storage of RE2 at each time period is either in zone 2 or zone 3, not in zone 1. Therefore, 
the model just assigns any value, which is either the lower bound or upper bound, as the 
solution of ai. 

The sensitivity analysis on the search ranges can provide water use and reservoir 
storage information for decision makers to decide new operating policies. The new rule 
curves and rationing factors should not only be able to improve the system performance, 
but also be easy to follow so that the reservoirs can be operated without difficulty. 
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For example, compare scenario 1 and 3, which have the same search range of 
rationing factors, the combination of new rule curves and rationing factors in scenario 3 
has a better objective function value, even though a: for RE1 of scenario 3 is smaller 
than that of scenarionl. It implies that the water use in scenario 3 may be smoother than 
scenario 1. 

The following is another example. Compare scenario 1 and 4, under the same rule 
curve search ranges, the objective function value of scenario 4 is better than that of 
scenario 1. For RE1 in scenario 4, it is easier to follow the hedging rules with the same ai 
and oi2 value. Note that since the performance of sensitivity analysis responds to the input 
data and initial conditions, it is difficult to establish relationships among these variables. 

In literature, it is generally difficult to deal with a MINLP problem. Fortunately, 
in this study, the transformation technique presented herein can take advantage of the 
special nonlinear formulation structure and convert the proposed MINLP model into a 
MILP model, which is easier to resolve by standard integer programming algorithms. 
However, if there is no such transformation technique existing, according to the analysis 
on other solution methodologies, the proposed hybrid SA may be suitable to handle the 
original MINLP model, even though the hybrid SA requires numerous evaluations and 
the computer time is much longer than for other methods. 



71 



Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 



Table 5-1 Upper Bounds and Lower Bounds Used When Searching for New Rule 
Curves for Simplified System 



Time 


RE1 


RE2 


Firm Curve 


Target Curve 


Firm Curve 


Target Curve 


LB 


UB 


LB 


UB 


LB 


UB 


LB 


UB 


1 


28 


32 


48 


50 


18 


22 


28 


30 


2 


25 


29 


45 


49 


15 


19 


25 


29 


3 


22 


26 


42 


46 


12 


16 


22 


26 


4 


19 


23 


39 


43 


9 


13 


19 


23 


5 


16 


20 


36 


40 


6 


10 


16 


20 


6 


13 


17 


33 


37 


3 


7 


13 


17 


7 


13 


17 


33 


37 


3 


7 


13 


17 


8 


16 


20 


36 


40 


6 


10 


16 


20 


9 


19 


23 


39 


43 


9 


13 


19 


23 


10 


22 


26 


42 


46 


12 


16 


22 


26 


11 


25 


29 


45 


49 


15 


19 


25 


29 


12 


28 


32 


48 


50 


18 


22 


28 


30 



Table 5-2 Upper Bounds and Lower Bounds Used When Searching for New Rationing 
Factors for Simplified System 



Rationing Factor 


RE1 


RE2 


LB 


UB 


LB 


UB 


cti 

0C2 


0.4 
0.6 


0.6 
0.8 


0.5 
0.7 


0.7 
0.9 
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Table 5-4 The Optimization Results of MINLP Model Solved by Transformation 
Technique for the Simplified System 



Time 


Water Supply 


Reservoir Beginning Storage 


DD1 


DD2 


RE1 


RE2 


DD2' 


DD2 2 


Storage 


Zone 


Storage 


Zone 


1 


8 


8 


9 


35 


2 


25 


2 


2 


8 


8 


9 


29 


2 


24 


2 


3 


8 


8 


10 


23 


2 


23 


3 


4 


6 


6 


9 


17 


1 


20 


2 


5 


8 


8 


10 


17 


2 


22 


3 


6 


8 


8 


10 


13 


2 


22 


3 


7 


6 


6 


10 


9 




22 


3 


8 


6 


6 


10 


11 




24 


3 


9 


6 


6 


10 


13 




26 


3 


10 


6 


6 


10 


15 




28 


3 


11 


6 


6 


10 


11 




28 


3 


12 


6 


6 


10 


7 




28 


3 



Table 5-5 The New Rule Curves for the Simplified System 



Time 


Ri 


El 


RE2 


Firm Curve 


Target Curve 


Firm Curve 


Target Curve 


1 


30 


50 


20 


30 


2 


29 


45 


15 


25 


-> 


23 


46 


12 


23 


4 


23 


43 


9 


21 


5 


17 


40 


6 


16 


6 


13 


37 


3 


13 


7 


13 


37 


7 


13 


8 


16 


40 


6 


16 


9 


19 


43 


9 


19 


10 


22 


42 


12 


22 


11 


25 


45 


15 


25 


12 


28 


48 


18 


28 
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Table 5-6 The New Hedging Rules for the Simplified System 



Reservoir Storage 



Rationing Factor 



RE1 



RE2 



Minimum Storage < S t < Firm Storage 

Firm Storage < S t < Target Storage 
Target Storage < S t < Maximum Storage 



60% 
80% 
100% 



70% 
90% 
100% 



S t : Beginning storage of RE1 (or RE2) in time period t. 



Table 5-7 Scenarios of Sensitivity Analysis on the Search Ranges of Rule Curves and 
Rationing Factors (Based on the Current Setting) 



Scenarios 


Rule Curves 


Rationin 


g Factors 


LB (Units) 


UB (Units) 


LB 


UB 


1 


-2 


+2 


-0.1 


+0.1 


2 


*> 

-j 
-4 


+3 
+4 


-0.1 
-0.1 


+0.1 
+0.1 


4 

5 


-2 
-2 


+2 
+2 


-0.2 
-0.3 


+0.2 
+0.3 



Table 5-8 Results of Sensitivity Analysis on Rationing Factors 



Scenarios 


R] 


El 


R] 


E2 


Objective 
Function 


ai 


a 2 


ai 


Ot2 


1 


0.60 


0.80 


0.70 


0.90 


-281.0 


2 
3 


0.60 
0.60 


0.77 
0.78 


0.50 
0.50 


0.80 
0.90 


-282.4 
-284.6 


4 
5 


0.70 
0.69 


0.70 
0.70 


0.80 
0.30 


1.00 
1.00 


-287.0 
-287.0 
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Figure 5- 1 Flowchart of v3-GBD 
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Figure 5-2 The Physical Scheme of Simulated Annealing 
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Figure 5-3 The Flowchart of Hybrid Simulated Annealing Approach 
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Figure 5-4 Evaluation of Objective Function of HSA1 and HSA2 
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Figure 5-5 The New Rule Curves for the Simplified System (Scenario 1 ) 
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Figure 5-6 The New Rule Curves for the Simplified System (Scenario 2) 
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Figure 5-7 The New Rule Curves for the Simplified System (Scenario 3) 
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Figure 5-8 The New Rule Curves for the Simplified System (Scenario 4) 
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Figure 5-9 The New Rule Curves for the Simplified System (Scenario 5) 
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CHAPTER 6 CASE STUDY: THE SOUTHERN 
REGIONAL SYSTEM, TAIWAN 



A real- world water distribution system, the southern regional water supply system 
in Taiwan, will be used as a case study to demonstrate the MILP model presented in 
Chapter 4 to optimize the system operation with reservoir rule curves and hedging rules. 
In addition, the MINLP model proposed in Chapter 5 is then used to evaluate and modify 
the current rule curves and hedging rules in this system. Figure 6-1 shows the 
geographical area of the region. Figure 6-2 shows layout of the water resources in this 
region. 

6.1 Description of the Southern Regional System 

The water distribution system of southern Taiwan covers the following two major 
river basins: the Tsengwen Basin and the Kaoping Basin. The system supplies water to 
the Chiayi, Tainan, Kaohsiung, and Pingtung Counties in southern Taiwan. A network 
representation of the distribution system is shown in Figure 6-3. The three major 
reservoirs in the Tsengwen Basin are the Tsengwen Reservoir, the Wusantou Reservoir, 
and the Nanhwa Reservoir. The Tsengwen Reservoir is the largest reservoir in Taiwan, 
and it is always operated jointly with the Wusantou Reservoir for flood control, and 
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agricultural, municipal and industrial water supply. Inflow to the Wusantou Reservoir is 
negligible. Hence, the two reservoirs can be lumped into a single reservoir called the T-W 
Reservoir. The Nanhwa Reservoir mainly provides municipal and industrial water supply. 
In the Kaoping Basin, several smaller reservoirs provide municipal and industrial water 
supply. The basic characteristics of all the reservoirs are shown in Table 6-1. In addition 
to the reservoirs mentioned, the system has many diversion structures. 

In the southern part of Taiwan, almost 90% of the annual rainfall occurs from 
May to October. As can be seen from Figure 6-4, the distribution of the average annual 
surface runoff is highly uneven. This creates a challenge to reservoir operators. For 
example, if a drought occurs and an insufficient amount of water is stored in the 
reservoirs, future water supply will be in jeopardy. To conserve water and minimize the 
impact of a severe water shortage in the future, the reservoirs must be operated under the 
hedging rules. 

The network representation in Figure 6-3 shows 14 source nodes, 27 diversion 
nodes. 21 demand nodes, five reservoir nodes, 13 junction nodes, one low flow 
requirement node, and 88 links. Reservoir nodes and diversion nodes are the two types of 
control nodes that control distribution of the water. A diversion node either transfers 
water to or from another basin, or it supplies demand. Undiverted water is discharged to 
the ocean. 

In Taiwan, for planning purposes, the municipal and industrial water demands are 
usually considered as one kind of demand, called the public demand. Therefore, two 
kinds of demands, the agricultural and public demands, are considered in this case study. 
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The T-W Reservoir has three rule curves, and the Nanhwa Reservoir has two rule 
curves. The rule curves and the hedging rules of the T-W Reservoir and the Nanhwa 
Reservoir are described below and shown in Tables 6-2-6-4 and Figure 6-5. 

The interpretation of the rule curves and of the hedging rules for the T-W 
Reservoir is as follows: 

1 . If the sum of the beginning storages of the Tsengwen and Wusantou Reservoirs is 
greater than that of the target curve, the planned demand at each demand node is met 
at 100%. 

2. If the sum of the beginning storages of the Tsengwen and Wusantou Reservoirs is 
greater than that of the firm curve and less than that of the target curve, the planned 
demand for the agricultural use is met at 75%; but the planned demand for public use 
is met at 100%. 

3. If the sum of the beginning storages of the Tsengwen and Wusantou Reservoirs is less 
than that of the firm curve, the planned demand for agricultural use is met at 50%; but 
the planned demand for public use is met at 80%. Note that the flood control curve 
applies only to the Tsengwen Reservoir. In our analysis, flood control reservation is 
strictly observed. 

The interpretation of the rule curves and of the hedging rules for the Nanhwa 
Reservoir is as follows: 

1 . If the beginning storage is greater than that of the target curve, the planned demand at 
each demand node is met at 100%. 
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2. If the beginning storage is greater than that of the firm curve and less than that of the 
target curve, the planned demand at each demand node is met at 90%. 

3. If the beginning storage is less than that of the firm curve, the planned demand for 
each demand node is met at 80%. 

Examining the network configuration reveals the following: demand nodes Dl 
and D4 can receive water from the Tsengwen Reservoir; nodes D2, D3, and D5 can 
receive water from the Wusantou Reservoir; and node D6 can receive water only from 
the Nanhwa Reservoir. Water supplies to these six demand nodes are required to satisfy 
the hedging rules in the optimization model. In the current setting, demand nodes D7 and 
D8 are not subject to the hedging rules. This is because D7 and D8 can receive water 
from all of the three reservoirs; and, additionally, two inflows (13 and 14) can also supply 
water to D7 and D8. 

Consistent with the irrigation practices in Taiwan, the time period for reservoir 
management and operation is 10 days. A month is divided into three time periods, in 
which the first two time periods always contain 10 days and the last time period varies 
from 8 to 1 1 days, depending upon the particular calendar month under consideration. 
The water year in Taiwan is from November 1 to October 3 1 . A water year invariably has 
36 time periods. In this study, the 10-day basis is preferred as the time period, 
recognizing that certain time periods may not be exactly 1 days. Additional assumptions 
made in the optimization model include the following: (1) the source of water supply 
considered is from the surface water; (2) no loss occurs during the flow transshipment; 
(3) evaporation loss from the reservoir is negligible; (4) return flow is negligible. 
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In the network representation, the one low flow requirement node, referred to as 
the low flow node, in the system can be handled by assigning a lower bound for the low 
flow arc (Figure 6-6(a)). During normal periods of inflow, this constraint is usually 
satisfied; however, during periods of drought, this constraint may cause infeasibility. To 
circumvent this, the low flow node is considered as a demand node in the optimization 
model, and a dummy node associated with the original low flow node is created (Figure 
6-6(b)). Accordingly, Figure 6-6(a) is transformed into Figure 6-6(b) with one dummy 
node and two artificial arcs. Note that the capacity of the arc, AL, in Figure 6-6(a) is also 
modified as shown in Figure 6-6(b). 

In the objective function, the low flow node is specified with a corresponding unit 
cost. Therefore, the optimization model will "supply" the low flow requirement as much 
as possible. During normal periods of inflow, there is sufficient water to satisfy the low 
flow requirement, and all surplus water will flow into the dummy node and then return to 
the low flow node. The total amount of water that flows into the low flow node is the 
summation of the low flow requirement and the surplus water. However, during a 
drought period, if there is a shortage of water to satisfy the low flow requirement, the 
amount of water that flows into the low flow node will be less than the low flow 
requirement, but infeasibility is avoided and no water will flow into the dummy node. 
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6.2 System Optimization with Current Rule Curves and Hedging Rules 

The following sections will describe the MILP model for the southern system in 
Taiwan in order to consider the reservoir rule curves and hedging rules and also to 
provide a discussion of the numerical results. 

6.2.1 MILP Optimization Model 

In the MILP model, the following operational objectives are considered herein: 
(1) maximization of water supply for public demand; (2) maximization of water supply 
for agricultural demand; (3) satisfaction of the low flow requirement; (4) maximization of 
reservoir storage. Using the weighting method combines the four objectives to reflect 
their relative importance. 

The most important objective is to distribute the available water and to maximize 
supply for public demand, or, in other words, to minimize any shortage in public demand. 
The second most important objective is to maximize the water supply for agricultural 
demand. The low flow requirement is the third concern. Least important, maximization of 
the reservoir storage, ensures that there is sufficient amount of water in storage for future 
water supply and continuing operations. Therefore, the objectives are prioritized in the 
following descending order: public demand, agricultural demand, low flow requirement, 
and reservoir storage. The maximization problem is converted to a standard minimum- 
cost problem by multiplying the objective function with a minus sign. 
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The one-year MILP model for the southern system in Taiwan can be formulated 
as follows: 

Min -t - KoU^- UfX*w- K c *),a*- 1(^x^1(6-1) 
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Eqs. (4-22)~(4-29) 

where 

Nfl = subset of N for reservoir nodes; 

N/ = subset of N for inflow nodes; 

Ng = subset of N for junction nodes; 

N/>v = subset of N for diversion nodes; 

N/, = subset of N for the low flow nodes; 

INj,, = amount of inflow to inflow node/; 

LWjj = low flow requirement at low flow node/; 
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Np = subset of N for the public demand nodes; 

N« = subset of N for the public demand nodes not subject to hedging rules; 

N.4 = subset of N for the agricultural demand nodes; 

No = subset of N for the agricultural demand nodes not subject to hedging rules; 

Cp= unit cost of supply to public demand; 

C A = unit cost of supply to agricultural demand; 

Cl = unit cost of flow to the low flow node; 

Cr= unit cost of storage arc of the reservoir. 

The composite objective function, Eq. (6-1), is a combination of four objectives. 
Eqs. (6-2)~(6-5) are the continuity equations for inflow, junction, diversion and reservoir 
nodes (here, b h , is zero). Eq. (6-6) specifies low flow requirement. Eq. (6-7) specifies the 
water supply constraints for demand nodes that are not subject to hedging rules. Note that 
the water supply is assumed to be equal to or less than the planned demand. The effect of 
this constraint is to avoid a situation under which certain demand nodes may receive a 
surplus supply and some other demand nodes may experience a severe shortage. Eqs. (4- 
22M4-29) have been defined before. Note that *=1, NK=\, P£*=l for D1-D6; K/={T-W 
Reservoir} for D1-D5, and K,={Nanhwa Reservoir} for D6. The weighting factors of 
these objectives are shown in Table 6-5. The values are obtained from consultations with 
engineers in Taiwan. In the minimum-cost optimization problem, the weighting factor for 
the public demand is the smallest to ensure that its priority is the highest. Similarly, the 
weighting factors of the other three objectives are assigned with values which correspond 
to their priorities. 
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6.2.2 Numerical Results and Discussions 

In this study, the MILP optimization problem is solved by LINGO. All the input 
data is provided by the Water Resources Bureau in Taiwan. The inflow data obtained for 
this study comprises 20 years of data (1974-1993). Figure 6-7 shows the variation of the 
annual inflows. Table 6-6 shows the annual inflow for each year. The third low inflow, 
which occurred in 1984 (5569.11 million cubic meters), is selected to test the 
optimization model. 

Note that the demand data for each demand node changes with the time period 
during the year. There are 36 different planned demands in 36 time periods for each 
demand node. The demand data for each demand node is shown in Table 6-7. The initial 
reservoir storage of each reservoir used in the optimization model is the average storage 
of the first 10-day period of the historical records. The initial storage of Tsengwen, 
Wusantou, and Nanhwa Reservoirs are 303.60, 45.07, and 1 10.65 million cubic meters, 
respectively. 

In the 36-time-period optimization model, there are 4,188 decision variables, 
including 210 integer variables and 3,137 constraints. Because the initial reservoir 
storage is known, it is not necessary to consider die integer variables for the beginning of 
the first time period. The input data consists of the inflow, planned demands, water 
supply limitation, link capacity limitation, and reservoir storage limitation for each time 
period. The reservoir rule curves and the hedging rules are also part of the input data. 
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The optimal regional water supply allocation solution has been obtained using 
LINGO. The results indicate that the reservoir release and water supply satisfy the rule 
curves and the hedging rules simultaneously. From Table 6-8, it is evident that the 
optimal solutions of the supply and reservoir storage zones for Dl to D5, the T-W 
Reservoir, D6, and the Nanhwa Reservoir satisfy the hedging rules. For example, the 
solution shows that the beginning storage of the T-W Reservoir in the second time step is 
in zone 2, and the corresponding optimal water supply for the Dl to D5 is in zone 2, not 
in zone 1 or zone 3. 

After testing the optimization model for the drought year, the inflow data of 20 
years, including the normal, wet, and dry years, is used to run the optimization model 
once per year to obtain the optimal solutions of water supply and reservoir storage. For 
planning purposes, it is not necessary to show the optimal water supply to each demand 
node and the optimal reservoir storage of each reservoir in each time period of each year. 
Instead, the average water supply and reservoir storage variations of the entire system for 
each of the 36 time periods in one year after the 20 runs are evaluated. The purpose of 
this procedure is to examine the average performance of the optimal operation of this 
system with hedging rules and various hydraulic situations. It is assumed that the demand 
of each demand node is the same in each year, but varies with the time period. 

Depending on the hydrologic data, the running time for each year varies from 3 
seconds to 4495 seconds on a PC with 1GHz Pentium III Processor and 1 GB RAM. The 
results of average supply and shortage from the 20-optimization runs are shown in 
Figures 6-8-6-10. In Figure 6-8, the axis represents the 20 different water years with 
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different inflows, and the ordinate in each figure represents different shortage 
information. The information presented allows decision makers to evaluate the average 
water supply and shortage conditions under a long-term operation when the system is 
operated with the hedging rules. For example, from Figure 6-10, it can be seen that the 
average shortage is more severe from the second to the eleventh time periods every year, 
which implies that the probability of water shortage during these periods may be higher 
than at other periods; and the water supply from other sources or a compromise between 
different demands may be needed in advance to avoid a shortage. 

Table 6-9 shows the sensitivity analysis results of total water supply and shortage 
to rationing factors. The values in the table are the averages of 20 runs. The relationship 
between the water supply and rationing factors can be explored via this table. For 
example, considering demand nodes D1-D4, the water supply of case 2-1 (rationing 
factors of the T-W Reservoir reduced by 0.1) is similar to that of case 3-1 (rationing 
factors of the T-W Reservoir increased by 0.1). As to D5, in which the planned demand is 
larger, the water supply of case 2-1 is less than that of case 3-1. Meanwhile, examining 
the results of D6, which is directly affected by the operation of the Nanhwa Reservoir, 
there is no major difference in water supply between case 2-1 and case 3-1. The results 
show that the rationing factors of the T-W Reservoir are not sensitive to D1~D4 and D6, 
but sensitive to D5. 

The following is another example. For D6, the water supply of case 2-2 (rationing 
factors of the Nanhwa Reservoir reduced by 0.1) is less than that of case 3-2 (rationing 
factors of the Nanhwa Reservoir increased by 0.1). Meanwhile, for Dl~ D4, the water 
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supply of case 2-2 is similar to that of case 3-2; however, for D5, the water supply of case 
2-2 is less than that of case 3-2. The results show that the rationing factors of the Nanhwa 
Reservoir are not sensitive to D1-D4, but sensitive to D5 and D6. Since sensitivity 
analysis is predicated on input data and a given set of base and initial conditions, it is 
difficult to establish relationships among the variables. But, it is interesting to note that 
most of the water supply results from case 1-1 to case 3-2 are superior to that of case 0, 
which uses the original set of rationing factors. This implies that, the current hedging 
rules should be modified to improve system operation. 

In this study, for demand nodes, D7 and D8, the sources are 13, 14, the T-W 
Reservoir, and the Nanhwa Reservoir. To maximize the water supply to D7 and D8, the 
optimization model will first allocate water from 13 and 14. If this amount of water is 
insufficient to meet the planned demands, the optimization model will then release the 
water from the reservoirs to supply D7 and D8 to minimize the objective function, 
because D7 and D8 are public demand nodes and have the highest priorities. 

It is important to note that after examining the water supply results, the demand 
nodes in this system can be approximately divided into two parts: D1-D8 which have 
more water shortage, and D9-D21 which have less water shortage. This is mainly 
because the inflow that can supply D9-D21 is more plentiful than that for D1-D8. 
Consequently, system managers may consider the feasibility of inter-basin water supply 
to supplement the water delivery to D1-D8. 
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63 Evaluation of Rule Curves and Hedging Rules 

In the southern system in Taiwan, the recent population and economic growth 
have caused increasing demands in water supply. In addition, the sediments dramatically 
have changed the volumes of the reservoirs since they were constructed. For instance, up 
to this point, the Wusantou Reservoir has already lost approximately 40% of its original 
volume. Therefore, to improve the performance of system operation, it is necessary to 
modify the current reservoir rule curves and hedging rules, which were determined 
beforehand in the planning stage and may not be suitable to the current situation. 

In this system, Tsengwen Reservoir is the only reservoir, which has the flood 
control function. The purpose of the reservoir flood control curve is to provide a storage 
space to store the water when a flood occurs and to reduce the downstream flooding 
damage. Analyzing the flood control function of a reservoir requires specific hydraulic 
aspects and a smaller time scale. Because this study focuses on the system optimization 
for planning purposes, modifying the flood control curve of the Tsengwen Reservoir is 
not included here. Consequently, there are four rule curves that need to be modified: the 
firm storage and target storage curves of T-W Reservoir, and the firm storage and target 
storage curves of Nanhwa Reservoir. In addition, the six rationing factors are modified as 
well in this case study. Based on the analysis on the MINLP model in Chapter 5, the 
transformation technique is chosen to solve the MINLP model for this system. 
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6.3.1 MINLP Optimization Model 

The formulation of the one-year MINLP model for the southern system in Taiwan 
is similar to the MILP model presented earlier, but the rule curves and rationing factors 
are unknown here. To compare the performances of new and current rule curves and 
hedging rules, the hydraulic characteristics of 1984 are chosen as an illustrative example. 

To avoid obtaining a new set of rule curves with randomly oscillatory 
distribution, a set of lower bounds and upper bounds for rules curves are assumed. The 
range assumption is also applied to search for new rationing factors which will eliminate 
any possibility of a large difference between the new ai and 0:2. 

For the T-W Reservoir, the lower bounds and the upper bounds of the firm 
storage and target storage are assumed to be the original values minus and plus 50 million 
cubic meters. The allowed minimum value for the firm storage curve is assumed as 20 
million cubic meters. If the assumed lower bound of target storage is less than the upper 
bound of firm storage, these two values are substituted by each other. The lower bounds 
and the upper bounds of rationing factors for the agricultural demands are assumed as the 
original values minus and plus 0.1. The rationing factors for the public demands are 
arbitrarily assigned without significant difference. The assumed lower bounds and upper 
bounds of rule curves and rationing factors are shown in Table 6-10 and 6-11. 

For the Nanhwa Reservoir, the lower bounds and the upper bounds of the firm 
storage and target storage are assumed to be the original values minus or plus 10 million 
cubic meters. The allowed minimum value for the firm storage curve is assumed as 2 
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million cubic meters. If the assumed lower bound of target storage is less than the upper 
bound of firm storage, these two values are substituted by each other. The lower bounds 
and the upper bounds of rationing factors are assumed to be the original values minus and 
plus 0.05. The assumed lower bounds and upper bounds of rule curves and rationing 
factors are shown in Table 6-10 and 6-11. 

Note that the initial storage of each reservoir is known beforehand. Unlike the 
MILP model above, the rule curves are unknown. Also, it is not certain which storage 
zone the initial reservoir storage belongs to in the beginning. Therefore, the integer 
variables for the initial stage may be considered in the MINLP model. In this case, after 
examining the initial storage, lower bound, and upper bound of rule curves in the 
beginning of the first time period for the Tsengwen, Wusantou, and Nanhwa Reservoirs, 
it is found that the integer variables for the T-W Reservoir need be considered. The 
integer variables for Nanhwa Reservoir in the beginning are ignored because the initial 
storage of Nanhwa is certainly in zone 2. 

6.3.2 Numerical Results and Discussions 

Using the transformation technique, there are 5912 continuous variables, 216 
integer variables, and 6378 constraints in the resulting MILP model. A new set of rule 
curves and hedging rules can be obtained by solving this model. The resulting objective 
function value is -1370634, which is superior to the original objective function value, - 
1331513, under current rule curves and hedging rules. The comparison of system 
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operation between the current and new rule curves and hedging rules are shown in Table 
6-12. This result demonstrates that the system can improve the performance of reservoir 
operation if the resulting new rule curves and hedging rules are incorporated in the 
system operation. 

In general, one reservoir has only one set of rule curves and hedging rules for 
long-term operation. However, this set of rule curves and hedging rules may not suit 
drought periods. Therefore, determining the rule curves and hedging rules for reservoir 
operation in dry years and normal years are of interest in this study. 

To determine the rule curves and hedging rules for dry years, it is necessary to 
select a set of hydraulic years to represent the situation in dry years. In this study, the five 
years with the lowest total annual inflows in the historical records are chosen to represent 
the dry years. These representative years are: 1980, 1993, 1984, 1991, and 1986. 

The one-year basis MINLP model is then implemented to obtain the new rule 
curves and hedging rules for each year. After five model runs, five sets of rule curves and 
hedging rules are achieved. In each set, there may be a number of curve oscillations 
which occur. The final rule curves and hedging rules for dry year operation are 
determined by averaging the value at each time period of each rule curve. Thus the curves 
can be smoother and easier to follow. The final rationing factors are similarly determined 
by averaging the five sets of rationing factors. The suggested rule curves and rationing 
factors for dry year operation are shown in Figure 6-1 1 and Tables 6-13 ~ 6-15. 

When the system is operated with the current rule curves and hedging rules, the 
optimization results show that the average annual water shortage in these five dry years is 
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697.76 million cubic meters. If the system is operated with the suggested rule curves and 
hedging rules, the optimization results show that the average annual water shortage in 
these five dry years is 650.58 million cubic meters. Examining these two situations, it is 
evident that the suggested rule curves and hedging rules for dry year operation can 
improve the system operation. 

Similarly, the sixth to fifteenth years in the twenty historical records sorted in 
ascending order are chosen to represent the normal years. These years are: 1988, 1987, 
1989, 1983, 1979, 1985, 1982, 1978, 1992, and 1976. After ten model runs, the final rule 
curves and rationing factors are determined by averaging these ten model results. The 
rule curves and rationing factors for normal year operation are shown in Figure 6-12 and 
Table 6-16 -6-18. 

When the system is operated with the current rule curves and hedging rules, the 
optimization results show that the average annual water shortage in these ten normal 
years is 434.90 million cubic meters. If the system is operated with the suggested rule 
curves and hedging rules, the optimization results show that the average annual water 
shortage in these ten normal years is 375.26 million cubic meters. Examining these two 
situations, it is evident that the suggested rule curves and hedging rules for normal year 
operation can improve the system operation. 

The result of a comparison between these two sets of rule curves and hedging 
rules for the T-W Reservoir shows that the target curve for dry year operation is generally 
higher than that for normal year operation; and the firm curve for dry year operation is 
generally lower than that for normal year operation. It implies that the model prefers to 
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let the hedging be active for dry year operation earlier than that for normal year 
operation, thus saving the water for future use. Because the water in a dry year is saved 
earlier than during a normal year, the firm curve of dry year operation can be lower than 
that of normal year operation. A similar situation also applies to the Nanhwa Reservoir 
but the target curve of dry year operation is closer to that of a normal year. 

The average annual water shortage under the current and suggested rule curves 
and hedging rules for dry and normal year operation are summarized in Table 6-19. 

In summary, the current rule curves and hedging rules can be modified via the 
proposed MINLP model to improve the current system operation. The results obtained 
can provide information which includes the reservoir operation and water usage in the 
system so that decision makers can decide optimal operating policies for planning 
purposes. 
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Table 6-1 Reservoir Characteristics, Southern Regional System, Taiwan 



Reservoir 


Active Storage 
(Million Cubic Meters) 


Main Purposes 


Tsengwen 

Wusantou 

Nanhwa 

Chengchinghwu 

Fengshan 


692.69 

83.76 

149.46 

5.00 

8.30 


F,A,I,M,P 

I,M 

A,I,M 

I,M 



F : Flood Control 
A : Agricultural Water Supply 
I : Industrial Water Supply 
M : Municipal Water Supply 
P : Hydropower 
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Table 6-2 Rule Curves of the T-W Reservoir and Nanhwa Reservoir 



Units: Million Cubic Meters 



Time 


T-W 


Nanhwa 


Firm 


Target 


Flood Control 


Firm 


Target 




Curve 


Curve 


Curve 


Curve 


Curve 


1 


240.0 


360.0 


500.0 


96.4 


122.4 


2 


230.0 


350.0 


480.0 


92.4 


111.4 


3 


220.0 


330.0 


460.0 


82.4 


104.8 


4 


210.0 


310.0 


440.0 


77.8 


96.4 


5 


200.0 


280.0 


420.0 


73.2 


88.4 


6 


190.0 


250.0 


400.0 


60.9 


84.4 


7 


170.0 


220.0 


380.0 


57.9 


80.4 


8 


160.0 


190.0 


360.0 


48.9 


71.4 


9 


150.0 


175.0 


340.0 


46.1 


66.0 


10 


120.0 


145.0 


320.0 


35.7 


54.9 


11 


90.0 


115.0 


300.0 


29.1 


47.4 


12 


80.0 


100.0 


280.0 


22.5 


40.9 


13 


55.0 


80.0 


260.0 


13.3 


32.2 


14 


40.0 


65.0 


240.0 


4.0 


23.4 


15 


30.0 


50.0 


220.0 


2.8 


18.0 


16 


30.0 


40.0 


220.0 


2.8 


27.8 


17 


30.0 


50.0 


220.0 


32.2 


96.4 


18 


40.0 


70.0 


230.0 


50.4 


149.5 


19 


60.0 


90.0 


250.0 


67.8 


149.5 


20 


80.0 


120.0 


300.0 


64.2 


149.5 


21 


105.0 


150.0 


360.0 


69.6 


149.5 


22 


130.0 


180.0 


420.0 


86.4 


149.5 


23 


155.0 


210.0 


460.0 


133.9 


149.5 


24 


180.0 


240.0 


581.0 


143.1 


149.5 


25 


200.0 


270.0 


617.0 


143.1 


149.5 


26 


240.0 


300.0 


617.0 


143.1 


149.5 


27 


280.0 


330.0 


617.0 


140.8 


149.5 


28 


280.0 


360.0 


617.0 


140.8 


149.5 


29 


280.0 


360.0 


617.0 


138.5 


149.5 


30 


280.0 


360.0 


617.0 


136.2 


149.5 


31 


280.0 


360.0 


617.0 


134.4 


143.8 


32 


280.0 


360.0 


580.0 


131.6 


142.2 


33 


280.0 


360.0 


570.0 


122.4 


140.6 


34 


280.0 


360.0 


560.0 


118.0 


131.6 


35 


260.0 


360.0 


540.0 


113.6 


129.3 


36 


250.0 


360.0 


520.0 


104.8 


127.0 
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Table 6-3 Hedging Rules of the T- W Reservoir 



Reservoir Storage 



Rationing Factor 



Agricultural Demand 



Public Demand 



Minimum Storage < S T , t + S w ,t < Finn Storage 

Firm Storage < Sr,t+ S w ,t < Target Storage 

Target Storage < S-r,t + Sw,t 



50% 
75% 
100% 



80% 

100% 

100% 



Si.t : Beginning storage of the Tsengwen Reservoir in time period t. 
Sw.t : Beginning storage of the Wusantou Reservoir in time period t. 



Table 6-4 Hedging Rules of the Nanhwa Reservoir 



Reservoir Storage 


Rationing Factor 


Minimum Storage < Sn.i < Firm Storage 

Firm Storage < S N ,t < Target Storage 

Target Storage < Sn.i 


80% 
90% 
100% 



Sn.i : Beginning storage of the Nanhwa Reservoir in time period t. 



Table 6-5 Weighting Factors Used for Southern Regional System, Taiwan 



Priority 



Objectives 



Order of 
Magnitude 



Weighting 
Factor 



1 
2 

3 
4 



Supply for the Public Demands 
Supply for the Agricultural Demands 

The Low Flow Requirement 
Reservoir Storage 



10 '-10' 
10" 2 ~10' 
10"'~10° 

iq-'-io 2 



1.0E+3 

1 
1.0E-3 
1.0E-6 
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Table 6-6 Annual Inflows, 1974-1993 (Sorted in Ascending Order) 



Year 



Annual Inflow (Million Cubic Meters) 



1980 


2840.80 


1993 


4082.37 


1984 


5569.11 


1991 


6525.27 


1986 


6877.09 


1988 


7068.27 


1987 


7228.66 


1989 


7444.10 


1983 


7454.69 


1979 


7873.59 


1985 


8189.84 


1982 


8591.37 


1978 


8615.88 


1992 


8984.99 


1976 


10015.27 


1974 


10562.96 


1981 


10599.54 


1975 


11051.27 


1990 


11583.94 


1977 


13316.01 
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Table 6-7 Demand Data 



Units: Cubic Meters Per Second 



Time 


Dl 


D2 


D3 


D4 


D5 


D6 


D7 


D8 


D9 


D10 


Dll 


1 


5.77 


2.93 


1.21 


0.01 


1.05 


9.45 


0.47 


2.12 


2.23 


2.33 


3.82 


2 


5.77 


2.93 


1.21 


0.02 


1.05 


9.45 


0.47 


2.12 


2.85 


3.22 


5.20 


3 


5.77 


2.93 


1.21 


0.05 


0.29 


9.45 


0.47 


2.12 


3.80 


3.87 


5.57 


4 


5.74 


3.11 


1.21 


0.07 


26.25 


9.31 


0.47 


2.56 


3.03 


3.35 


4.63 


5 


5.74 


3.05 


1.55 


0.08 


22.74 


9.31 


0.47 


2.56 


3.70 


4.70 


4.03 


6 


5.74 


3.05 


1.55 


0.12 


7.24 


9.51 


0.48 


2.12 


3.67 


3.87 


3.17 


7 


5.68 


3.05 


1.36 


0.10 


2.16 


9.51 


0.48 


2.12 


4.01 


5.11 


3.67 


8 


5.68 


2.93 


1.36 


0.09 


8.82 


9.51 


0.48 


2.12 


4.13 


4.73 


4.04 


9 


5.68 


2.93 


1.36 


0.08 


18.64 


9.51 


0.48 


2.12 


4.80 


5.12 


7.95 


10 


5.88 


2.93 


1.36 


0.09 


17.15 


9.51 


0.48 


2.12 


6.36 


6.46 


6.07 


11 


5.88 


2.93 


1.36 


0.05 


12.74 


9.45 


0.47 


2.12 


4.43 


4.33 


6.23 


12 


5.88 


2.93 


1.36 


0.01 


1.05 


9.45 


0.47 


2.58 


1.43 


1.43 


5.59 


13 


5.75 


3.05 


1.36 


0.00 


16.84 


9.45 


0.47 


2.58 


1.52 


1.52 


8.37 


14 


5.75 


3.05 


1.21 


0.00 


17.55 


9.45 


0.47 


2.58 


0.30 


1.30 


7.77 


15 


5.75 


3.05 


1.21 


0.00 


20.90 


9.51 


0.48 


2.58 


0.00 


0.00 


3.39 


16 


5.78 


3.05 


1.21 


0.00 


6.41 


9.31 


0.47 


2.58 


0.00 


0.00 


5.50 


17 


5.78 


2.98 


1.21 


0.00 


10.28 


9.31 


0.47 


2.58 


0.10 


0.00 


8.86 


18 


5.78 


2.98 


1.21 


0.00 


21.58 


9.31 


0.47 


2.18 


0.70 


0.95 


6.90 


19 


5.80 


2.98 


1.21 


0.01 


27.69 


9.31 


0.47 


2.18 


0.71 


0.75 


6.11 


20 


5.80 


3.05 


1.21 


0.00 


39.13 


9.31 


0.47 


2.33 


0.96 


1.56 


7.11 


21 


5.80 


3.05 


1.21 


0.01 


61.68 


9.45 


0.47 


2.33 


3.74 


3.99 


7.84 


22 


5.65 


3.05 


1.36 


0.00 


26.19 


9.45 


0.47 


2.33 


0.42 


0.42 


6.61 


23 


5.65 


3.11 


1.36 


0.00 


41.92 


9.45 


0.47 


2.33 


0.39 


0.39 


7.08 


24 


5.65 


3.11 


1.36 


0.00 


46.02 


9.45 


0.47 


2.58 


0.53 


0.53 


9.50 


25 


5.47 


3.11 


1.36 


0.00 


41.09 


9.41 


0.47 


2.58 


0.72 


0.72 


7.90 


26 


5.47 


3.11 


1.36 


0.02 


26.91 


9.41 


0.47 


2.12 


1.02 


1.52 


7.38 


27 


5.47 


3.05 


1.36 


0.00 


28.92 


9.41 


0.47 


2.12 


1.41 


1.67 


5.42 


28 


5.27 


3.05 


1.36 


0.05 


34.63 


9.41 


0.47 


2.12 


1.70 


1.79 


0.05 


29 


5.27 


3.05 


1.36 


0.06 


17.82 


9.41 


0.47 


2.12 


1.94 


2.94 


0.35 


30 


5.27 


2.98 


1.36 


0.04 


9.48 


9.51 


0.48 


2.12 


0.59 


0.59 


0.31 


31 


5.55 


2.98 


1.36 


0.06 


0.97 


9.51 


0.48 


2.12 


1.51 


1.79 


0.39 


32 


5.55 


2.98 


1.36 


0.07 


0.97 


9.51 


0.48 


2.12 


1.68 


1.68 


0.52 


33 


5.55 


3.05 


1.36 


0.08 


0.97 


9.51 


0.48 


2.12 


1.73 


1.95 


0.14 


34 


5.47 


3.05 


1.36 


0.10 


0.97 


9.51 


0.48 


2.12 


1.63 


1.63 


0.60 


35 


5.47 


3.05 


1.36 


0.09 


0.97 


9.51 


0.48 


2.12 


1.80 


1.77 


0.84 


36 


5.47 


3.05 


1.36 


0.08 


1.06 


9.51 


0.48 


2.12 


1.52 


1.72 


0.60 | 
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Table 6-7 Demand Data (Continued) 

















Units: Cubic Meters Per Seconc 


Time 


D12 


D13 


D14 


D15 


D16 


D17 


D18 


D19 


D20 


D21 


1 


1.43 


3.94 


3.40 


9.08 


3.57 


2.80 


0.54 


2.12 


2.89 


7.08 


2 


1.43 


3.94 


8.19 


9.08 


3.57 


4.01 


0.54 


3.07 


2.89 


9.04 


3 


1.43 


4.33 


11.14 


9.08 


3.57 


3.46 


0.54 


3.69 


2.89 


11.34 


4 


1.43 


3.94 


8.83 


9.34 


3.57 


3.48 


0.57 


3.62 


2.89 


9.20 


5 


1.43 


3.94 


8.83 


9.34 


3.57 


3.50 


0.57 


4.50 


3.25 


9.51 


6 


1.43 


3.15 


9.56 


9.34 


3.57 


3.50 


0.57 


5.29 


3.25 


10.51 


7 


1.43 


3.94 


8.75 


9.31 


3.57 


3.51 


0.54 


5.18 


3.25 


9.90 


8 


1.43 


3.94 


8.75 


9.31 


3.57 


3.47 


0.54 


5.40 


3.25 


10.21 


9 


1.43 


4.33 


8.67 


9.31 


3.57 


3.47 


0.54 


5.38 


3.25 


10.17 


10 


1.43 


3.94 


8.88 


9.64 


3.57 


3.24 


0.43 


7.01 


2.89 


11.83 


11 


1.43 


3.94 


8.88 


9.64 


3.57 


3.24 


0.43 


6.37 


2.89 


6.20 


12 


1.43 


3.94 


7.20 


9.64 


3.57 


1.90 


0.43 


1.96 


2.89 


4.83 


13 


1.43 


3.94 


3.55 


9.38 


3.57 


0.42 


0.49 


2.12 


2.89 


5.99 


14 


1.43 


3.94 


0.94 


9.38 


3.57 


0.06 


0.49 


0.48 


2.89 


1.25 


15 


1.43 


4.33 


0.00 


9.38 


3.57 


0.00 


0.49 


0.00 


3.25 


0.59 


16 


1.43 


3.94 


0.00 


9.65 


3.57 


0.00 


0.44 


0.00 


3.25 


0.86 


17 


1.43 


3.94 


2.88 


9.65 


3.57 


0.46 


0.44 


0.55 


3.25 


2.44 


18 


1.43 


3.94 


4.17 


9.65 


3.57 


2.95 


0.44 


1.98 


2.89 


5.12 


19 


1.43 


3.94 


7.59 


10.13 


3.57 


3.23 


0.54 


3.62 


2.89 


7.61 


20 


1.43 


3.94 


9.73 


10.13 


3.57 


2.58 


0.54 


3.18 


2.89 


8.06 


21 


1.43 


4.33 


11.51 


10.13 


3.57 


2.63 


0.54 


4.16 


2.89 


10.93 


22 


1.43 


3.94 


11.04 


10.27 


3.57 


2.43 


0.45 


2.33 


3.25 


5.68 


23 


1.43 


3.94 


11.04 


10.27 


3.57 


2.43 


0.45 


2.17 


3.25 


4.95 


24 


1.43 


4.33 


10.66 


10.27 


3.57 


2.43 


0.45 


2.82 


3.25 


6.45 


25 


1.43 


3.94 


11.60 


10.37 


3.57 


2.63 


0.53 


1.22 


3.25 


6.28 


26 


1.43 


3.94 


11.60 


10.37 


3.57 


2.63 


0.53 


0.00 


3.25 


1.77 


27 


1.43 


3.94 


11.03 


10.37 


3.57 


1.70 


0.53 


0.00 


2.89 


0.75 


28 


1.43 


3.94 


7.88 


10.06 


3.57 


0.35 


0.51 


0.03 


2.89 


0.96 


29 


1.43 


3.94 


3.38 


10.06 


3.57 


0.46 


0.51 


0.20 


2.89 


1.37 


30 


1.43 


4.33 


0.00 


10.06 


3.57 


0.00 


0.51 


0.06 


2.89 


0.75 


31 


1.43 


3.94 


0.00 


9.77 


3.57 


0.24 


0.55 


0.19 


2.89 


0.92 


32 


1.43 


3.94 


0.00 


9.77 


3.57 


0.42 


0.55 


0.59 


2.89 


1.74 


33 


1.43 


3.94 


0.00 


9.77 


3.57 


0.77 


0.55 


1.21 


2.89 


2.97 


34 


1.43 


3.94 


0.00 


9.47 


3.57 


0.23 


0.42 


1.61 


2.89 


4.74 


35 


1.43 


3.94 


0.00 


9.47 


3.57 


0.30 


0.42 


1.77 


2.89 


5.10 


36 


1.43 


4.33 


0.00 


9.47 


3.57 


0.95 


0.42 


1.02 


2.89 


3.94 
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Table 6-8 The Optimal Zones of the Water Supply and Reservoir Beginning Storage 
(1984) 



Time 


T-W 


Nanhwa 


Dl 


D2 


D3 


D4 


D5 


Reservoir 


D6 


Reservoir 


1 


2 


2 


2 


2 


2 


2 


2 


2 


2 


2 


2 


2 


2 


2 


2 


2 


2 


3 


2 


2 


2 


2 


2 


2 


2 


2 


4 


2 


2 


2 


2 


2 


2 


2 


2 


5 














2 


2 


6 














2 


2 


7 














2 


2 


8 














2 


2 


9 














2 


2 


10 


2 


2 


2 


2 


2 


2 


2 


2 


11 


2 


2 


2 


2 


2 


2 


2 


2 


12 


1 


1 


1 


1 


1 


1 


2 


2 


13 


3 


3 


3 


3 




3 


2 


2 


14 


3 


3 


3 


3 


3 


3 


3 


3 


15 


3 


3 


3 


3 


3 


3 


3 


3 


16 


2 


2 


2 


2 


2 


2 


3 


3 


17 


2 


2 


2 


2 


2 


2 


2 


2 


18 


2 


2 


2 


2 


2 


2 


2 


2 


19 


2 


2 


2 


2 


2 


2 


2 


2 


20 


3 


3 


3 


3 


3 


3 


3 


3 


21 


2 


2 


2 


2 


2 


2 


3 


3 


22 














3 


3 


23 














3 


3 


24 














3 


•^ 
j 


25 














3 


3 


26 














3 


3 


27 














3 


3 


28 














3 


3 


29 














3 


3 


30 














2 


2 


31 














3 


3 


32 














2 


2 


33 














2 


2 


34 














3 


3 


35 














2 


2 


36 














2 


2 1 



109 



Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 



en 

O 
CO 

u. 

00 

c 

"c 
o 

« 
0* 



u 
ao 

co 

r 
o 



a 

3 
C/3- 
>_ 

S5 

o 

c/i 
'33 
_>> 

<n 

c 

< 



c/s 

C 

u 
c/3 



OS 

I 

es 



la 






j 














"~r » m o 


— r- 


00 — 




s 




*!--> ® 


T m 


fl o 




u 


vO 

n 




! 












! 








X! 

3 




3 
t*) 


VOjin — 


S 2 


r- ^r 




u 

c 
o 




'TiO ci 
vo:Ov os 
fN j(N CN 


Ov 00 
00 00 
IN (N 


O r»i 
Os OS 
cm rvi 




1 

in 


Q 


en 


00 
3 


<N 00 

ri *r 

"£ 2 


IN CI 
O — 
in Ov 
SO CN 


OS ~ 




c 

3 


3 

in 


00 
OS 


t— — 
CN o 


t- so 
f- VO 


VO fN 
Ov OS 








r- 
ci 


ci — 

©N CI 
C") «r 


O vd 
r- o 


o' r-^ 
t — 








Si 
1/3 


— ir- — 


vo ooii^i r— 

CN CNjtN CN 








©j© o 


©' ojo © 










i 


1 
1 










I 








3 


m 00 T 
Ov|© O 


ON r-|tN 00 
Ov ON © OS 








Ol 


© Oj— ©' 

i 








C/3 


ooioo •v 
— ivO vo 


Ov CI 1 r— Ov 

v© voir- r- 








— ;© O 


o ojo o" 






Q 
























3 


vOjvO O 
©im »o 


in — ir- m 
m r-it *r 








C/3 


o'io' © 


o o'io' ©' 








C/i 


voir- m 

00 ! OS Ov 


ci r— i ci on 

O void CN 








— io o 


— Ol 






«N 

Q 




: : i 


















3 


ci|cn tt 

— io o 


vo cNjr- o 

©s cijvo r- 








</) 


CMIci r*i 
Ov : Os Ov 


cn ciitN rsi 

Ov Os 1 Ov Ov 








-C 
C/5 


00 
00 


vo ci 

O OS 


cn vo 


r- in 






__ 


fi 


«N — 


CN — 


<N <N 






Q 


3 


CN 

Ov 


tt r- 
r- oe 


00 T 
v© r*i 


iri m 








r- 


rn ci 

r- r- 


r- r- 


r- r- 






u 

IO 

a 
U 


Ol ■ ■ ] 


— <N 









JZ 
C/3 


00 
vO 


Ov 

o 


o 
o 


o 
o 


o 
o 


O O 

o o 




o 


o 


o 


o 


o 


o o 


00 

Q 
































3 


00 
OS 


r- 


VO 
VO 


VO 
vo 


vO 
vO 


VO vo 
VO vo 




(/5 


Ov 
VO 


o 
r- 


o 
r- 


o 
r- 


o 
r- 


O ©' 

r- r- 






00 


o 
o 


o 
o 


o 
o 


o 
o 


o o 
o o 




o 


o 


o 


o 


o 


o o 


r- 
Q 
































3 


00 


vO 

r- 


VO 

r- 


vO 

r- 


SO 

f- 


VO vo 

r- r- 




V) 


T 


Tf 


■*r 


rr 


T 


rr ■«• 


u 
a 

U 


O 


- 




i 


(N 

i 

tN 


— CM 

■ i 
m ci 



>. 








o. 




o. 




-t 


<u 


CA 


Oil 




CO 


u 


C 


(0 


o 

x: 






to 


to 






o 


o 


H 


r- 



e/3 t/3 



to 

Si 

c 

^^ CQ 03 CS CQ CQ CQ 

<- * S * i i * 

£££££££ 
w"" 03 w w n Co CQ 

2 Z 2 Z Z Z 2 

Sr u. h. Um u w b. 

|U i2 a a a a & 

~ 00 O Ov 00 0v O 

°^ o" — " o d d — 

<=> I! II II II II II 

If r* <N r*t m r* r>* 

g a a 8 a a a 

^ t~ os" oo* r-^ oo" oC 
>" d d d d d d 

9f ll_ ll_ II II II ll_ 

*~. 3, SI 3, & 3, 3, 

1 jt 5?" ^* 3=" £ 5:" 

El i till 
u H H f- f- E- r- 

"° -a" -a ■a* -a" -o" t»* 

^Zm h (9 (0 q (B a 

^EEEEEE 

2 u u u u u u 

O. 13 -a -ra -o -a -a 

o .y .y .y .2 .y .y 

^Z Z Z I Z ": 

3 3 3 3 3 3 
_^ o. q. o. a a. a. 

"T" i- If I- k. U L. 
II OOOOOO 

" IB lis IS IS IS IS 

a '"v /•■ *> ' — v /— ^ /-s /^ 

. Os O OS p O O 

d ii^ T n m 7 IT T 
^ a a a a a a 

"-^ 1— os r— 00 os 00 

-LT o" d o" d d ©' 

St ||_ ||_ ||_ n_ 11 11 

1 3=" ^* 3:' 3=* £ 5t' 

Ei 1 1 l 1 t 
„ H H H H t- (- 

"° T3* -a -o* -o" -a" ■a* 

-= c c c c c c 

2 ts to 3 « to eg 
= E E E E E E 

£Z CU Q> O 4) u u 

3 T3 -O -O T3 -O ■O 

'5, S "S 2 S 2 S 

"3 3 3 3 3 3 

^ .y .y .y .- .y .s 

in eb ob eb eb ab eb 

(— - to tQ to CO CQ CS 
,_ : h- k- t— u u h. 
^ ,0 ^ ^ .0 ,0 

|| IW t»« Im Ita lM Ita 

8 in in m in m in 
- vo 00 vo r- 00 r- 

ij- d o" dodo 

<~> 11 11 11 11 11 11 

1 r* n r* n n *s 

- a a a a a a 

^ ^* vo* ^* m* vo" m* 

O d o" d © ©' p' 

c 11 u_ 11 11 11 T 

£j a a a a a a 

3 '*-' ^^ ^^ N-' ^-^ N-' 

ii — cn ~ rsi — rsi 

1 1 1 1 1 1 

O — — CN (N r»i r*> 

ft> U tU 4> U U U 

3 2 en to wj to to 

cs co co ca cd cb 

U U U U O U U 



110 



Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 



Table 6-10 Upper Bounds and Lower Bounds Used When Searching for New Rule 
Curves for T-W Reservoir and Nanhwa Reservoir 



Time 


T-W Reservoir 


Nanhwa Reservoir 


Firm Curve 


Target Curve 


Firm Curve 


Target 


Curve 


LB 


UB 


LB 


UB 


LB 


UB 


LB 


UB 


1 


190.0 


290.0 


310.0 


410.0 


86.4 


106.4 


112.4 


132.4 


2 


180.0 


280.0 


300.0 


400.0 


82.4 


101.4 


102.4 


121.4 


3 


170.0 


270.0 


280.0 


380.0 


72.4 


92.4 


94.8 


114.8 


4 


160.0 


250.0 


260.0 


360.0 


67.8 


86.4 


87.8 


106.4 


5 


150.0 


230.0 


250.0 


330.0 


63.2 


78.4 


83.2 


98.4 


6 


140.0 


200.0 


240.0 


300.0 


50.9 


70.9 


74.4 


94.4 


7 


120.0 


170.0 


220.0 


270.0 


47.9 


67.9 


70.4 


90.4 


8 


110.0 


140.0 


210.0 


240.0 


38.9 


58.9 


61.4 


81.4 


9 


100.0 


125.0 


200.0 


225.0 


36.1 


55.1 


56.1 


76.0 


10 


70.0 


95.0 


170.0 


195.0 


25.7 


44.7 


45.7 


64.9 


11 


40.0 


65.0 


140.0 


165.0 


19.1 


37.4 


39.1 


57.4 


12 


30.0 


50.0 


130.0 


150.0 


12.5 


30.9 


32.5 


50.9 


13 


20.0 


30.0 


105.0 


130.0 


3.3 


22.2 


23.3 


42.2 


14 


20.0 


30.0 


90.0 


115.0 


2.0 


13.0 


14.4 


33.4 


15 


20.0 


30.0 


80.0 


100.0 


2.0 


10.8 


12.8 


28.0 


16 


20.0 


30.0 


80.0 


90.0 


2.0 


12.8 


17.8 


37.8 


17 


20.0 


30.0 


80.0 


100.0 


22.2 


42.2 


86.4 


106.4 


18 


20.0 


30.0 


90.0 


120.0 


40.4 


60.4 


139.5 


149.5 


19 


20.0 


40.0 


110.0 


140.0 


57.8 


77.8 


139.5 


149.5 


20 


30.0 


70.0 


130.0 


170.0 


54.2 


74.2 


139.5 


149.5 


21 


55.0 


100.0 


155.0 


200.0 


59.6 


79.6 


139.5 


149.5 


22 


80.0 


130.0 


180.0 


230.0 


76.4 


96.4 


139.5 


149.5 


23 


105.0 


160.0 


205.0 


260.0 


123.9 


139.5 


143.9 


149.5 


24 


130.0 


190.0 


230.0 


290.0 


133.1 


139.5 


149.5 


149.5 


25 


150.0 


220.0 


250.0 


320.0 


133.1 


139.5 


149.5 


149.5 


26 


190.0 


250.0 


290.0 


350.0 


133.1 


139.5 


149.5 


149.5 


27 


230.0 


280.0 


330.0 


380.0 


130.8 


139.5 


149.5 


149.5 


28 


230.0 


310.0 


330.0 


410.0 


130.8 


139.5 


149.5 


149.5 


29 


230.0 


310.0 


330.0 


410.0 


128.5 


139.5 


148.5 


149.5 


30 


230.0 


310.0 


330.0 


410.0 


126.2 


139.5 


146.2 


149.5 


31 


230.0 


310.0 


330.0 


410.0 


124.4 


133.8 


144.4 


149.5 


32 


230.0 


310.0 


330.0 


410.0 


121.6 


132.2 


141.6 


149.5 


33 


230.0 


310.0 


330.0 


410.0 


112.4 


130.6 


132.4 


149.5 


34 


230.0 


310.0 


330.0 


410.0 


108.0 


121.6 


128.0 


141.6 


35 


210.0 


300.0 


310.0 


410.0 


103.6 


119.3 


123.6 


139.3 


36 


200.0 


300.0 


310.0 


410.0 


94.8 


114.8 


117.0 


137.0 
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Table 6-1 1 Upper Bounds and Lower Bounds Used When Searching for New Rationing 
Factors for T-W Reservoir and Nanhwa Reservoir 



Rationing Factors 


T-W Reservoir 


Nanhwa Reservoir 


Agricultural Demand 


Public Demand 


LB 


UB 


LB 


UB 


LB 


UB 


a-> 


0.40 
0.65 


0.60 
0.85 


0.75 
0.85 


0.85 
1.00 


0.75 
0.85 


0.85 
0.95 



Table 6- 1 2 Comparison of System Operation Between the Current and New Rule 
Curves and Hedging Rules (1984) 



Rule Curves 

and Hedging 

Rules 


Objective 

Function 

Value 


Number of 

10-day 

Shortage 


Water 
Supply 
(MCM) 


Water 
Shortage 
(MCM) 


Percentage 

of Shortage 

(%) 


Current 
New 


-1331513 
-1370634 


36 
29 


2061.95 
2236.31 


630.13 

455.77 


23.41 
16.93 
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Table 6-13 The Suggested Rule Curves for Dry Year Operation 



Units: Million Cubic Meters 



Time 


T-W 


Nanhwa 


Firm 


Target 


Flood Control 


Firm 


Target 




Curve 


Curve 


Curve 


Curve 


Curve 


1 


190.00 


390.00 


500.00 


86.40 


112.40 


2 


200.00 


320.35 


480.00 


90.00 


103.25 


3 


210.00 


300.00 


460.00 


84.40 


96.29 


4 


178.00 


272.20 


440.00 


75.24 


87.80 


5 


182.00 


266.00 


420.00 


69.28 


83.20 


6 


152.00 


252.00 


400.00 


58.90 


74.40 


7 


130.00 


230.00 


380.00 


55.90 


73.44 


S 


116.00 


210.00 


360.00 


46.90 


61.40 


9 


105.00 


205.00 


340.00 


46.30 


56.10 


10 


80.00 


175.00 


320.00 


33.30 


45.70 


11 


45.00 


140.00 


300.00 


26.42 


42.76 


12 


38.00 


134.00 


280.00 


19.86 


32.50 


13 


22.00 


110.00 


260.00 


18.42 


23.30 


14 


24.00 


90.00 


240.00 


6.40 


14.40 


15 


22.00 


81.68 


220.00 


5.52 


12.80 


16 


26.00 


81.95 


220.00 


6.32 


21.80 


17 


24.00 


84.00 


220.00 


34.20 


86.40 


18 


22.00 


102.00 


230.00 


44.40 


141.50 


19 


24.00 


110.00 


250.00 


62.80 


139.50 


20 


38.00 


130.00 


300.00 


62.20 


139.50 


21 


64.00 


191.00 


360.00 


79.60 


139.50 


22 


80.00 


220.00 


420.00 


84.40 


139.50 


23 


105.00 


249.00 


460.00 


130.14 


143.90 


24 


130.00 


254.00 


581.00 


135.66 


149.50 


25 


164.00 


292.00 


617.00 


136.94 


149.50 


26 


202.00 


314.00 


617.00 


136.94 


149.50 


27 


240.00 


336.28 


617.00 


130.80 


149.50 


28 


238.20 


339.38 


617.00 


130.80 


149.50 


29 


232.23 


335.54 


617.00 


132.90 


148.50 


30 


240.29 


338.51 


617.00 


126.20 


146.86 


31 


230.00 


362.00 


617.00 


128.16 


144.40 


32 


230.00 


362.00 


580.00 


132.20 


141.60 


33 


230.00 


362.00 


570.00 


119.68 


132.40 


34 


246.00 


362.00 


560.00 


113.44 


128.00 


35 


210.00 


350.00 


540.00 


114.30 


123.60 


36 


220.00 


350.00 


520.00 


102.80 


117.00 
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Table 6-14 The Suggested Hedging Rules of the T-W Reservoir for Dry Year Operation 



Reservoir Storage 


Rationing Factor 


Agricultural 
Demand 


Public 
Demand 


Minimum Storage < Sf.t + Sw.t < Firm Storage 

Firm Storage < S T , t + Sw.t < Target Storage 

Target Storage < Sr,t + Sw,t 


44% 
66% 
100% 


83% 
98% 
100% 



S-p.t : Beginning storage of the Tsengwen Reservoir in time period t. 
Sw.t : Beginning storage of the Wusantou Reservoir in time period t. 



Table 6-15 The Suggested Hedging Rules of the Nanhwa Reservoir for Dry Year 
Operation 



Reservoir Storage 



Rationing Factor 



Minimum Storage < S N , t < Firm Storage 

Firm Storage < Sn.i < Target Storage 
Target Storage < Sw.t 



85% 
95% 
100% 



Sn.i : Beginning storage of the Nanhwa Reservoir in time period t. 
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Table 6-16 The Suggested Rule Curves for Normal Year Operation 



Units: Million Cubic Meters 



Time 


T-W 


Nanhwa 


Firm 


Target 


Flood Control 


Firm 


Target 




Curve 


Curve 


Curve 


Curve 


Curve 


1 


210.00 


390.00 


500.00 


101.40 


112.40 


2 


240.00 


310.00 


480.00 


99.50 


102.40 


3 


240.00 


280.00 


460.00 


90.40 


94.80 


4 


223.00 


260.00 


440.00 


84.54 


87.80 


5 


206.00 


250.00 


420.00 


76.88 


83.20 


6 


182.00 


240.00 


400.00 


68.90 


74.40 


7 


150.00 


220.00 


380.00 


65.90 


72.40 


8 


128.00 


210.00 


360.00 


56.90 


61.40 


9 


115.00 


200.00 


340.00 


55.10 


56.10 


10 


85.00 


170.00 


320.00 


42.80 


45.70 


11 


55.00 


140.00 


300.00 


35.57 


39.10 


12 


42.00 


130.00 


280.00 


29.06 


32.50 


13 


26.00 


105.00 


260.00 


22.20 


23.30 


14 


26.00 


90.00 


240.00 


11.90 


14.40 


15 


26.00 


80.00 


220.00 


9.92 


12.80 


16 


26.00 


80.00 


220.00 


11.72 


17.80 


17 


26.00 


80.00 


220.00 


40.20 


86.40 


18 


26.00 


90.00 


230.00 


58.40 


139.50 


19 


32.00 


113.00 


250.00 


75.80 


139.50 


20 


54.00 


134.00 


300.00 


72.20 


139.50 


21 


55.00 


200.00 


360.00 


79.60 


139.50 


22 


110.00 


190.00 


420.00 


94.40 


139.50 


23 


138.00 


216.00 


460.00 


139.50 


143.90 


24 


149.22 


242.00 


581.00 


138.22 


149.50 


25 


192.00 


257.00 


617.00 


138.22 


149.50 


26 


226.00 


296.00 


617.00 


139.50 


149.50 


27 


237.23 


331.03 


617.00 


130.80 


149.50 


28 


232.06 


330.00 


617.00 


130.80 


149.50 


29 


231.07 


330.00 


617.00 


129.60 


148.70 


30 


230.44 


330.00 


617.00 


128.86 


146.72 


31 


263.00 


346.00 


617.00 


132.86 


144.40 


32 


254.00 


338.00 


580.00 


132.20 


141.60 


33 


246.00 


338.00 


570.00 


128.78 


132.40 


34 


238.00 


338.00 


560.00 


120.24 


128.00 


35 


228.00 


320.00 


540.00 


119.30 


123.60 


36 


210.22 


320.00 


520.00 


112.80 


117.00 
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Table 6-1 7 The Suggested Hedging Rules of the T-W Reservoir for Normal Year 
Operation 



Reservoir Storage 


Rationing Factor 


Agricultural 
Demand 


Public 
Demand 


Minimum Storage < Sr,t+ Sw.t < Firm Storage 

Firm Storage < S T , t + S w ,t < Target Storage 

Target Storage < Sr.t + Sw.t 


45% 
69% 
100% 


84% 
98% 
100% 



Sj.t : Beginning storage of the Tsengwen Reservoir in time period t. 
Sw.t : Beginning storage of the Wusantou Reservoir in time period t. 



Table 6- 1 8 The Suggested Hedging Rules of the Nanhwa Reservoir for Normal Year 
Operation 



Reservoir Storage 



Rationing Factor 



Minimum Storage < Sn.i < Firm Storage 

Firm Storage < S N ,t < Target Storage 
Target Storage < Sn.i 



81% 
95% 
100% 



Sn.i : Beginning storage of the Nanhwa Reservoir in time period t. 



Table 6-19 Average Annual Water Shortage Comparison of Current and Suggested Rule 
Curves and Hedging rules 



Rule Curves and 
Hedging Rules 



Average Annual Water Shortage (MCM) 



Dry Year 



Normal Year 



Current Setting 
Suggested Setting 



697.76 
650.58 



434.90 
375.26 
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Figure 6-1 Area Covered by the Southern Regional System, Taiwan 



Taiwan Strait 
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Figure 6-2 Layout of Water Resources in the Southern Regional System, Taiwan 
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Figure 6-4 Average Annual Inflow in the Southern Region 
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Figure 6-5 Rule Curves of the T-W Reservoir and Nanhwa Reservoir 



(a) The Rule Curves of T-W Reservoir 
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(b) The Rule Curves of Nanhwa Reservoir 
-.*..- Firm Storage - -• - - Target Storage 
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Figure 6-6 Network Transformation for Low Flow Node 
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Figure 6-7 Annual Inflow in the Southern Region 
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Figure 6-8 The Optimized Results 



(a) Number of 10-day Shortage 




(b) Percentage of Shortage 
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Figure 6-9 The Optimized Long-term Reservoir Storage 
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Figure 6- 1 The Optimized Long-term Demand and Supply 
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Figure 6-1 1 The Suggested Rule Curves for Dry Year Operation 



(a) Rule Curves of T-W Reservoir 



• Flood Control 
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Figure 6-12 The Suggested Rule Curves for Normal Year Operation 
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CHAPTER 7 CONCLUSIONS AND FUTURE RESEARCH 



Reservoir management and operation play an important rule in overall water 
resources planning and management. During normal periods of operation, allocating the 
available water from different reservoirs to meet the planned demands imposed by 
competing users can be achieved without difficulty. However, during periods of drought, 
or when anticipating a drought, water demand cannot be met 100% and water shortages 
occur. The reservoir rule curves and the hedging rules provide operating guidelines that 
can be followed to minimize the impacts of a drought when water in the reservoirs are 
depleted because of insufficient inflows. 

This study focuses on incorporating the rule curves and hedging rules in the 
optimization of reservoir management and operation. Furthermore, compared with the 
current setting, a superior set of rule curves and hedging rules are achieved to improve 
the current system operation. The conclusions and future research of this study will be 
described in the following sections. 

7.1 Conclusions 

Listed below are a number of conclusions reached by this study: 
1 . A multi-reservoir, multi-purpose MILP model that considers both the rule curves and 
the hedging rules is developed for optimizing a large-scale regional water distribution 
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system. The developed model was tested and verified on a simplified system. Then 
this model was successfully applied to the regional water distribution system in the 
southern region of Taiwan. The results obtained demonstrated the applicability and 
utility of the proposed MILP model. A manuscript summarizing the results has been 
accepted for publication (Tu et al., 2002). 

2. A MINLP model is developed by expanding the proposed MILP model to achieve the 
new reservoir rule curves and hedging rules to improve the current reservoir 
operation. By applying this model to the simplified system, the new rule curves and 
hedging rules obtained are demonstrated to be superior. Then this methodology was 
successfully applied to the regional water distribution system in the southern region 
of Taiwan to achieve new rule curves and hedging rules. The results showed that the 
performance of system operation could be improved if the new rule curves and 
hedging rules are implemented. In addition, two types of rule curves and hedging 
rules, dry year and normal year, are proposed for realistic multi-reservoir system 
operation for different hydrology situations. 

3. To solve the presented MINLP model which is a nonconvex optimization problem, a 
number of methodologies are proposed. These methodologies are: Generalized 
Benders Decomposition, Penalty Coefficient Methods, Pseudo-integer Method, 
Simulated Annealing, and one mathematical transformation technique. The 
advantages and drawbacks of these methodologies are pointed out in the previous 
discussion. A comparison of these methodologies reveals that the transformation 
technique takes advantage of the mathematical structure and is most suitable to 
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solving the presented MINLP model. 

4. In order to determine a set of rule curves and hedging rules for a multi-purpose 
reservoir, numerous affecting factors such as the system hydrology characteristics, 
reservoir operating purposes, water supply for planned demands, and orderly rule 
curves, need to be considered. Consequently, not only a mathematical optimization 
model is necessary for decision making, but also the establishment of a cooperative 
relationship between the system managers, reservoir operators and water users. 

5. When a set of rule curves and hedging rules are determined, they may not be 
applicable for the reservoir after years of operation. Therefore, it is necessary to 
modify the rule curves and hedging rules from time to time to optimize the reservoir 
operation. 

7.2 Future Research 

Several works for future research are suggested based on the discussion of this 
study: 

1 . Future studies may consider the hydraulic uncertainty to analyze the reliability of new 
rule curves and hedging rules for reservoir management and operation. 

2. The use of ground water may be considered to supplement the surface water supply. 

3. More system characteristics may be taken into account to improve the proposed 
optimization models. These characteristics include: reservoir evaporation, reservoir 
storage change with time, loss/gains during water transshipment, and so on. 
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4. In general, the pattern of hedging consists of two or three rules curves with 
corresponding rationing factors for a multi-purpose reservoir. This kind of pattern 
may need to be reevaluated (such as introducing four rule curves) for better reservoir 
operation performance. 
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